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1.  OBJECTIVE  OF  THE  PROGRAM 


This  Report  summarizes  the  two  year  research  program  focussed  on  the  gradual 
loss  of  the  mechanical  strength  and,  eventually,  integrity  of  concrete  structures 
exposed  to  chemically  aggressive  substances  commonly  found  in  the  environment. 
This  program  is  a  continuation  of  the  research  program  on  the  micro-mechanics  of 
concrete  conducted  previously  by  the  principal  investigator  and  his  associates. 
Both  of  these  programs  were  sponsored  by  the  Civil  Engineering  Program, 
Directorate  of  Aerospace  Sciences  of  the  United  States  Air  Force  Office  of 
Scientific  Research. 

The  objective  of  these  studies  was  primarily  aimed  on  the  better  under¬ 
standing  of  the  complex  physico-chemical  phenomena  leading  to  the  deterioration 
of  concrete  reacting  with  the  water-borne  chemicals.  The  research  described  in 
this  Report  is  obviously  of  basic  nature  since  it  represents  an  inquiry  into  the 
basic  phenomena.  Yet  it  also  has  undeniably  practical  implications  as  well.  It 
appears  obvious  that  this  type  of  research  can  later  serve  as  a  basis  for  the 
improvements  of  the  current  design  practices  ir.  order  to  enhance  the  economy  and 
reliability  of  concrete  structures  such  as  runway  pavements. 

The  nature  of  this  problem  and  its  scope,  spanning  over  several  scientific 
disciplines,  requires  parallel  investigations  of  distinctly  different,  but  inherently 
coupled,  processes.  Some  of  this  processes  are  seldom,  if  ever,  considered 
within  the  context  of  mechanics  of  concrete  structures.  In  an  effort  to  make  this 
Report  self-sufficient  and,  consequently,  more  useful  for  the  reader  the  analytical 
modelling  of  chemical  reactions  and  diffusion  processes  is  in  the  sequel  replicated 
in  a  sufficient  detail.  A  list  of  references  containing  an  even  more  comprehensive 
inquire  into  these  processes  was  provided  to  help  the  interested  reader  feeling  the 
need  for  a  better  understanding  of  the  phenomena. 

In  discussing  a  sequence  of  phenomena  as  complex  as  the  one  reported 
herein  it  was  often  necessary  to  focus  on  a  specific  problem.  A  relative  loss  of 
generality  was  more  than  compensated  by  the  depth  of  the  insight  gained  through 
considerations  of  the  selected  case.  The  relative  paucity  of  specific  experimental 
data  mandated  a  not  entirely  supported  selection  of  a  most  likely  scenario 
describing  the  process. 
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The  general  nature  of  the  proposed  model  is  micromechanical.  The  ultimate 
benefit  arising  from  such  an  approach  is  that  it  enables  establishment  of  rational 
relationships  between  the  chemical  composition  of  the  hardened  cement  paste, 
morphology  of  the  pore  system,  and  defect  distribution  on  one  hand  and  the 
macroscopic  (phenomenological)  attributes  of  the  specimens  such  as  strength  and 
durability  on  the  other  hand. 

The  considered  problem  of  chemo-micro-mechanics  of  concrete  is,  obviou¬ 
sly,  too  wide  to  be  addressed  in  a  comprehensive  manner  within  a  two  year 
research  program.  Even  though  most  of  the  considerations  to  be  discussed  in  the 
sequel  are  of  general  validity  (or  will,  at  the  very  least,  become  a  part  of  a  more 
general  analytical  model)  the  studies  summarized  in  this  Report  are  focussed  on: 

-  a  specific  do»ble-decomposition  chemical  process  leading  to  the  formation 
of  ettringite  resulting  from  the  exposure  of  the  hardened  mortar  to  water¬ 
borne  magnesium  sulphates, 

-  isothermal  processes  in  order  to  reduce  the  number  of  thermodynamic  var¬ 
iables  and  thermal  diffusion, 

-  small  to  moderate  pore  pressures  in  order  to  be  able  to  neglect  the  barro- 
diffusion  and  the  influence  of  the  pore  pressure  on  the  stress-strain  relation¬ 
ship  in  the  solid  skeleton  of  the  hardened  cement  paste,  and 

-  the  pre-critical  processes  with  regard  to  the  diffusion  and  mechanical  res¬ 
ponse. 

A  generalization  of  these  studies  to  the  critical  and  post-critical  phenomena, 
pressure  driven  flow  through  a  porous  media,  and  many  other  interesting 
phenomena  is  expected  to  become  a  part  of  subsequent  studies. 


2.  INTRODUCTION 


The  general  objective  of  this  report  is  to  relate  the  chemical  composition  of  the- 
concrete  and  the  morphology  of  its  microstructure  to  its  chemo-mechanical 
response  when  exposed  to  chemically  aggressive  ambients.  It  is  in  this  sense  that 
the  label  chemo-micro-mechanics  seems  both  descriptive  and  justified. 

Naturally,  this  work  is  by  no  means  the  first  one  with  such  a  lofty 
objective.  For  example,  the  chemically  assisted  fracture  was  studied  by  Rice 
(1978),  Charles  (1978),  Cherepanov  (1374),  Fuller  and  Thomson  (1979), 
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Thomson  (1980),  and  many  others  on  different  levels  of  sophistication.  A  rather 
general  discussion  of  degradation  and  chemomechanics  was  published  by  Aifantis 
(1980,  1981)  who  should  also  be  credited  with  an  initial  continuum  study  of 
embrittlement  (Colios  and  Aifantis,  1982).  Many  other  papers  with  similar 
objective  can  be  located  in  the  existing  literature.  It  is,  nevertheless,  felt  that  this 
Report  is  not  only  different  from  the  similar  efforts  but  that  it  also  establishes  a 
comprehensive  framework  for  studies  of  an  entire  class  of  problems.  Finally,  this 
work  has  a  practical  aspect  as  well  since  it  aspires  to  become  ultimately  a  basis 
for  a  powerful  tool  in  design  of  concrete  structures. 

Concrete  is  a  complex  composite  material  containing  several  active 
ingredients  capable  of  reacting  with  the  aggressive  chemicals  commonly  found  in 
the  environment.  These  reactions  are  often  deleterious  for  the  integrity  of 
concrete  causing  distress  of  structures  and  necessitating  costly  repairs  and  even 
costlier  replacements. 

As  pointed  out  by  Moavenzadeh  (1971)  the  failure  of  a  structure  is  not  in 
itself  a  physical  phenomenon  but  a  condition  or  state  of  a  complex  process 
characterized  by  reduced  reliability  and  inability  of  the  structure  to  perform  its 
functions.  Consequently,  in  order  to  be  able  to  predict  the  failure  it  is  necessary 
to  analyze  the  processes  leading  to  the  distress  and  establish  criteria  needed  to 
distinguish  the  conditions  portending  the  onset  of  failure.  Thus,  the  analytical 
model  must  describe  the  evolution  of  variables  representing  physically 
identifiable  measures  of  state  and  formulate  conditions  defining  the  limiting 
response  of  the  material.  For  example,  brittle  failure  of  concrete  occurs  either  as 
a  result  of  a  runaway  growth  of  a  single  crack  of  preferential  geometry  or  due  to 
coalescence  (cooperative  action)  of  microcracks  into  a  macrocrack.  Hence, 
appropriate  internal  variables  (Krajcinovic  1989)  must  be  introduced  into  the 
analytical  model  to  describe  the  manner  in  which  the  microcracks  affect  the 
mechanical  response  of  the  solid  (compliance,  diffusivity,  etc.).  For  low  to 
moderate  crack  densities  a  single  internal  variable  (volume  average)  suffices. 
However,  at  the  localization  threshold  the  direct  interaction  becomes  a  dominant 
attribute  and  an  appropriate  internal  length  variable  must  be  introduced  in  order  to 
be  able  to  predict  the  critical  phenomena.  At  the  threshold  of  the  critical  regime 
(distress)  the  evolution  of  the  variable  describing  the  state  of  damage  (local  loss 
of  integrity)  is  rapidly  accelerated.  From  the  physical  standpoint  the  process 
becomes  non-local.  Matnemahcally,  the  problem  oecomes  strongly  nonlinear  and 


in  the  case  of  the  localization  the  governing  differential  equations  cease  to  be 
elliptical  causing  further  difficulties  in  analyses. 

Analytical  modeling  of  the  environmental  degradation  of  the  hardened 
cement  paste  involves  considerations  of  a  complex,  and  not  totally  understood, 
physico-chemical  process  interrelating  diffusion,  adsorption,  a  double 
decomposition  chemical  reaction  and  evolution  of  damage  (micro  and  macro¬ 
cracking).  The  complexity  of  these  processes  can  be  overwhelming  since  they  in 
an  essential  manner  depend  on  the  chemical  composition  and  microstructure  of 
concrete.  It  is,  therefore,  not  surprising  that  a  reasonably  comprehensive,  and 
convincing,  analytical  modelling  of  this  phenomenon  is  not  as  yet  available. 
Moreover,  the  experimental  evidence  is  often  sketchy  and  directed  towards 
conventional  phenomenological  modelling  being  of  limited  utility  in 
micromechanical  studies. 

Among  many  interesting  physico-chemical  processes  this  paper  focuses  on 
the  external  sulphate  attack  on  concrete.  Contained  in  a  wide  range  of  industrial 
and  agricultural  effluents  the  magnesium  sulphate  proves  to  be  a  potent  cause  of 
degradation  of  a  wide  range  of  concrete  structures  such  as  pavements,  culverts, 
bridges,  hydraulic  structures,  etc.  The  internal  sulphate  attack  (attributable  to 
sulphates  within  the  aggregates,  see  Ouyang,  et  al.  1988)  will  not  be  considered 
within  this  Report.  Naturally,  a  reasonably  straightforward  modification  of  the 
proposed  model  should  suffice  to  enable  analyses  of  the  internal  attack  as  well. 

The  degradation  of  concrete  as  a  result  of  exposure  has  been  recognized  by 
the  ACI  Committee  201  which  in  its  "Guide  to  Durable  Concrete"  -  ACI  201. 2R- 
77  reflects  on  this  problem  in  a  rather  cursory  manner.  In  particular,  the 
Committee  suggests  that  in  structures  exposed  to  sulphate  attack  (seawater  and 
groundwaters  in  western  and  norther  U.S.)  only  the  type  V  cement  (less  than  5% 
of  C3A  or  20%  of  C4AF  +  2C3A)  be  used.  It  also  recommends  classification  of 
exposure  into  four  classes  with  respect  to  the  sulphate  concentration.  The 
Committee  further  recommends  lower  w/c  ratios  of  0.45  and  advises  against  use 
of  calcium  chlorides  (reducing  sulphate  resistance) 

The  sulphates  can  damage  concrete  also  by  a  purely  physical  action  by 
crystallization  in  the  pores  of  concrete.  The  growing  crystals  (not  reacting  with 
the  cement)  can  exert  sufficient  pressures  to  cause  spalling  of  the  outer  layers  of 
concrete  (Reading  1975).  Naturally,  even  the  sulohate-resistant  Type  V  cements 
are  not  impervious  to  this  type  of  attack. 


The  objective  of  this  study  is  to  formulate  a  micro-mechanical  theory 
modeling  the  dominant  aspects  of  the  underlying  processes,  blurring 
simultaneously  the  unnecessary  details.  The  most  important  aspect  of  this  study 
consists  in  providing  a  rational  explanation  of  the  experimentally  measured 
trends.  Specifically,  the  primary  concern  centers  on  providing  a  rational 
analytical  model  relating  chemical  composition  of  concrete  and  its  morphology 
with  the  mechanical  response.  For  example,  it  has  been  clearly  established  that 
the  volumetric  expansion  of  concrete  depends  not  only  on  the  concentration  of 
sulphates  and  the  time  of  exposure,  but  also  on  the  composition  of  concrete  and 
in  particular  on  the  mass  fraction  of  the  tricalcium  aluminate.  Even  though  the 
basic  problem  is  very  complex  whenever  possible  and  advisable  an  effort  will  be 
made  to  achieve  the  objective  and  retain  the  simplicity  and  tractability  without 
compromising  the  rigor  and  physical  insight  into  the  underlying  phenomena. 


3.  ANALYTICAL  MODEL 


Consider  a  hardened  cement  paste  specimen  saturated  with  clean,  chemically 
inert,  water,  which  is  at  time  t  =  0  immersed  into  a  pool  of  water  containing  a 
known,  and  constant,  concentration  of  magnesium  sulphate  MgS04-  The  solute 
diffuses  into  the  hardened  cement  paste  specimen,  triggers  a  sequence  of  complex 
chemical  reactions  and  forms  an  expansive  reaction  product  (ettringite).  The 
attendant  local  stress  concentrations  in  the  material  surrounding  the  inclusion  in 
which  the  reaction  took  place,  are  often  sufficient  to  cause  microcracking  of  the 
hardened  cement  paste.  The  basic  thermodynamic  variables  needed  for  the 
description  of  the  process  during  which  the  mechanical  strength  of  the  concrete 
exposed  to  the  magnesium  sulphates  (or  some  similar  aggressive  chemical 
substances)  degrades,  are: 

-  solute  concentration  c(x,t)  of  the  magnesium  sulphate  in  the  diffusing 
water, 

-  concentration  of  the  adsorbent  (gypsum)  cg(x,t)  on  the  adsorbate  (wall  of 
the  pores  of  hardened  cement  paste), 

-  concentration  of  the  reaction  product  (ettringite)  ce(x,t)  formed  during  the 
double  decomposition  reaction, 

-  elastic  strain  e(x,t)  in  the  specimen,  and 


-  mechanical  history  lecording  parameter  (microcrack  density)  co(x,t). 

The  rates  of  change  of  these  thermodynamic  variables  are  related  through 
the  kinetic  equations  such  as  those  governing  the  rate  of  chemical  reactions, 
diffusion  equation,  adsorption  isotherms  and  damage  evolution  laws.  The 
coefficients  in  these  equations  must  be  determined  through  careful  experiments  in 
conjunction  with  micromechanical  considerations  neither  of  which  are,  according 
to  the  consulted  literature,  commonly  available. 

Finally,  it  is  necessary  to  determine  the  thermodynamic  state  of  the  mate¬ 
rial  as  defined  by  the  two  state  parameters:  diffusivity  and  mechanical 
compliance.  In  both  cases  the  consideration  was  at  this  point  restricted  to  the  pre- 
critical  regime.  Nevertheless,  the  derivation  of  micromechanically  rigorous 
expressions  for  the  diffusivity  and  compliance  tensors  requires  consideration  of 
the  internal  damage,  i.e.  presence  of  a  large  number  of  microcracks  distributed 
over  a  large  part  of  the  affected  volume.  In  a  subsequent  research  program,  to  be 
proposed  shortly,  these  investigations  will  be  extended  to  the  post-critical 
regimes  considering  the  propagation  of  the  percolation  and  fragmentation  fronts 
through  the  hardened  cement  paste  and  the  behavior  of  material  swept  by  these 
two  fronts. 

The  processes  of  diffusion,  adsorption,  chemical  reaction  and  micro- 
cracking  are  inherently  coupled.  For  example,  the  extent  of  the  chemical  reaction 
obviously  depends  on  the  transport  of  magnesium  sulphates  through  the  hardened 
cement  paste  (diffusion  or  flow  through  the  porous  media).  The  rate  of  transport 
of  aggressive  chemical  substances  is  enhanced  by  microcracking  which  opens 
new,  wider  and  more  direct  paths  for  the  percolating  fluid.  However,  the 
microcrack  nucleation  and  evolution  is  directly  attributable  to  the  formation  of 
ettringite  crystals  which  are  in  turn  the  products  of  the  chemical  reactions. 

Many  other  coupling  effects  such  as  influence  of  the  macro-stress  sign 
(compression  or  tension)  on  the  crack  width  and,  therefore,  diffusivity,  influ¬ 
ence  of  the  pore  pressure  on  the  stress-strain  relationship,  thermal  effects 
associated  with  the  exothermic  reactions,  etc.  will  at  this  point  be  considered  as 
being  of  secondary  importance  and,  consequently,  ignored. 

In  summary,  the  analytical  model  will  consist  of  several  partial  and 
ordinary  differential  equations  defining  the  salient  aspects  of  the  problem  such  as 
the  rate  of  chemical  reactions,  diffusion  with  adsorption  and  the  mechanical 
equilibrium.  These  equations  will  be  formulated  in  this  Report  sequentially  and 
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almost  independent  of  each  other.  The  Report  will  also  contain  a  description  of 
simplifications  introduced  to  enhance  the  tractability  and  brief  discussions  of 
trends  in  future  development  of  the  proposed  model.  Naturally,  the  accuracy  and 
predictive  ability  of  the  model  will  be  tested  against  the  available  experimental 
data.  The  main  emphasis  of  this  comparison  will  be  in  relating  the  macro¬ 
response  of  the  mortar  specimen  to  its  micro-structure  and  chemical  composition. 

The  flow  chart  shown  in  Fig.  1  represents  the  general  structure  of  the 
model  to  be  developed  in  detail  in  the  subsequent  sections  of  this  study. 


Fig.  1.  Schematic  structure  of  the  proposed  model  for  the  sulphate  corrosion 

induced  degradation  of  concrete 
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4.  CHEMICAL  REACTION 


A  review  of  the  existing  literature  reveals  a  reasonable  concensus  with  regard  to 
the  main  points  of  the  so-called  double  decomposition  reaction  leading  to  the 
formation  of  ettringite  during  the  exposure  of  the  hardened  cement  paste  to  the 
water  borne  sulphates.  However,  the  exact  details  of  these  reactions  are  still  a 
matter  of  debate  and  disagreement.  In  this  paper  the  so-called  Hansen' s  model 
(1968)  will  be  taken  as  the  starting  point.  The  proposed  model  can  readily 
accommodate  other  decomposition  reactions  which  can  be  found  in  the  literature, 
without  changes  in  its  basic  structure. 

Recently  reported  experiments  (Odler  and  Gassser  1988,  Moukwa  1990) 
confirm  the  essential  nature  of  the  micromechanical  changes  associated  with 
MgSC>4  attack  on  the  cement  paste.  Careful  measurements  indicate  a  substantial 
reduction  of  the  total  porosity  (10%  -  20%,  Moukwa  1990)  resulting  from  the 
expansion  products  filling  the  pores.  The  expansion  develops  rather  slowly  since 
the  pore  space  suffices  to  accommodate  the  ettringite  crystals  initially.  After 
protracted  exposure  (in  excess  of  half  a  year)  the  internal  fracture  (damage) 
becomes  apparent  since  the  remaining  pore  space  is  insufficient  to  provide  enough 
room  for  the  expanding  needle-like  ettringite  crystals.  According  to  Moukwa 
(1990)  the  "examination  of  fractured  surfaces  in  the  scanning  electron  microscope 
indicated  that  the  formation  of  such  reaction  products  probably  was  the  main 
cause  of  deterioration." 

According  to  Hansen  (1968)  the  double  decomposition  reaction  taking  place 
in  the  hardened  cement  paste  inundated  by  a  solution  of  magnesium  sulphate  can 
be  summarized  by  the  following  equations 


Ca(OH)2 

+ 

MgS04  + 

2^0 

— > 

CaS04-2H20 

+ 

Mg(OH) 

m :  74.0 

+ 

120.4  + 

36 

172.1 

+ 

58.3 

v :  33.2 

+ 

72.1  + 

36 

— > 

74.1 

+ 

24.5 

and 


m  : 

v  : 


3Ca0Al203  +  3(CaS04-2H20)  +  26H20 


270.2 

88.8 


+ 


3(172.1)  +  26(18) 
3(74.1)  +  26(18) 


3Ca0Al203-3CaS04-32H20  (4.2) 
1254.5 
725.1 


+ 


where  the  molar  masses  m  (g/mol)  and  molar  volumes  v  (cm3)  are  listed  below 
each  chemical  substance  taking  part  in'  the  reaction  (see  also  Appendix  A). 

The  first  of  the  two  reactions  (4.1)  is  a  through-solution  reaction.  The 
reaction  products  (gypsum  and  magnesium  hydroxide)  are  readily  accommodated 
by  the  porosity  of  the  mortar.  Simple  computations  clearly  indicate  "that  concrete 
cannot  be  caused  to  expand  and  crack  by  the  simple  mechanism  of  filling  pores 
with  solids  by  a  through-solution  process"  (Hansen  1968, page  32). 

The  second  reaction  (4.2),  involving  tricalcium  aluminate,  newly  formed 
gypsum  and  water  needed  to  form  the  ettringite,  is  often  assumed  to  be  of  the 
solid-liquid  or  topochemical  kind.  In  this  case,  "the  space  available  locally  at  the 
surface  is  not  big  enough  to  accommodate  the  newly  formed  crystal"  (Soroka, 
1979,  page  153).  The  expansion  of  the  newly  formed  ettringite  crystal  may, 
therefore,  cause  large  local  tensile  stresses  sufficient  to  induce  cracking  of  the 
hardened  cement  paste  in  the  vicinity  of  the  pores. 

The  Hansen's  model  (4.1  -  2)  is  not  the  only  one  to  be  found  in  the  litera¬ 
ture.  Mather  (1968)  proposes  two  alternatives:  a  single  decomposition  and  a 
double  decomposition  (somewhat  similar  to  Hansen's)  reaction.  Biczok  (1972) 
proposes  a  change  in  the  reaction  (4.2)  and  recognizes  two  alternatives  of  that 
reaction  for  pH  below  or  above  12.5.  Mindess  (1981)  also  agrees  with  the  first  of 
the  two  reactions  but  allows  for  transformation  of  the  ettringite  to  monosulphate 
in  the  case  when  all  sulphates  are  consumed  before  the  tricalcium-aluminates  are 
completely  hydrated.  According  to  Lea  (1970)  and  Soroka  (1979)  the  second  of 
the  two  reactions  involves  C4AH19 . 

Naturally,  it  is  very  difficult  to  determine  which  one  of  these  reactions  is 
typical  of  the  considered  case.  Moreover,  on  the  technical  scale  the  ASTM  C1012 
Standard  recommends  the  aqueous  magnesium  sulphate  (MgS04-7H20)  to  be  used 
in  the  expansion  tests  (instead  of  MgS04)  The  listed  ASTM  Code  also  specifies 
the  range  of  the  pH  factor  between  six  and  eight. 


4.1.  Stoichiometric  analysis 

The  masses  of  the  formed  gypsum  and  ettringite  and  the  masses  of  the 
substances  "consumed"  in  the  course  of  the  chemical  reactions  (4. 1)  and  (4.2)  can 
be  determined  using  conventional  stoichiometric  analyses. 

Denote  by: 


11 


c(x,t)  -> 

concentration  of  magnesium  sulphate  in  the  water 

<p 

accessible  porosity  of  the  hardened  cement  paste 

m 

— > 

molar  mass  (g/mol) 

M 

— > 

total  mass  of  a  substance  in  volume  V  =  1  m3. 

V 

— » 

molar  volume  (cm3). 

y 

— * 

mas;:  densities 

ise  the  follow 

ing  subscripts; 

S' 

for 

magnesium  sulphate 

cm 

for 

cement 

c 

for 

calcium  hydroxide 

a 

for 

tricalcium  aluminate 

g 

for 

gypsum 

e 

for 

ettringite. 

Superscripts  c  and  /  stand  for  masses  consumed  and  formed  in  the 
course  of  a  given  chemical  reaction. The  relevant  data  for  the  molar  masses  and 
volumes  of  all  substances  involved  in  reactions  (4. 1 )  and  (4.2)  are  assembled  and 
tabulated  in  the  Appendix  A  to  this  Report. 

Consider  now  the  chemical  reaction  (4. 1).  The  mass  of  the  gypsum  formed 
during  this  reaction  may  be  controlled  either  by  the  supply  of  sulphates  (if  the 
ratio  Ms/ Mc  <  ms/mc  =  1.63)  or  by  the  available  mass  of  calcium  hydroxide  (if 
Ms/Mc>  1.63).  Rather  elementary  stoichiometric  analyses  suffice  to  determine  the 
mass  of  the  newly  formed  gypsum 

=  min  [(mg  /  ms)<f>  cysoiV  ,  (mg  /  mc)  fc  Mcm]  (4j) 

while  the  mass  of  the  calcium  hydroxide  consumed  in  the  course  of  the  reaction 
(4. 1 )  is 


Mcc  =  min  [(mc  /  ms)<p  cysol  V  ,  fc  Mcm  ]  (44) 

In  the  two  above  equations  fc  denotes  the  mass  fraction  of  calcium  hydroxide 
with  respect  to  the  mass  of  cement,  Ms  =  c<PysoiV  is  the  mass  of  sulphates  with 
Ysot  denoting  the  density  of  magnesium  sulphate  solution.  The  mass  fraction  of  the 


calcium  hydroxide  fc  c  an  be  readily  determined  by  a  similar  stoichiometric 
analysis  o.  the  hydration  reaction  of  the  C3S  and  C2S  .  Since  the  experiments 
indicate  no  dependence  of  the  volumetric  expansion  on  the  available  mass  of  the 
calcium  hydroxide  this  altern?.  ve  can  be,  at  the  moment,  safely  discarded. 

The  mass  of  the  consumed  tricalcium  aluminate  during  the  second  (topo- 
chemical)  reaction  (4.2)  can  be  then  determined  to  be 

Ma  =  mi '  [(ma  /  3 mg)  Mfg  ;  fa  Mcm  ]  (4.5) 

depending  again  on  whether  the  reaction  is  controlled  by  the  supply  of  gypsum 
(Mg/Ma<  1.91)  formed  during  the  preceding  through-solution  reaction,  or  the 
available  mass  of  the  tricalcium  aluminate  (when  the  inequality  sign  is  reversed). 
In  (4.5)  fa  denotes  the  mass  fraction  of  tricalcium  aluminate  calculated  with 
respect  to  the  mass  of  cement. 

In  conjunction  with  (4.3)  the  expression  (4.5)  can  be  rewritten  in  the  form 

Ma  =  min  [(ma  /  3ms)<p  c  YsoiV  ;  ( ma/3mc)f'cM  cm  \  fa  Mcm]  (4.5)! 

Finally,  the  mass  of  the  ettringite  formed  in  the  course  of  the  topochemical  re¬ 
action  (4.2)  is 

=  (me/  ma)  Ma  =  min  [(me  /  3ms)<p  cysolV  ;  (me  /  3mc)fc  Mcm  ;  (me  /  ma)fa  Mcm] 

(4.6) 

Consequently,  the  mass  of  the  ettringite  formed  is  controlled  by:  (a)  the  mass  of 
the  available  sulphates,  (b)  available  calcium  hydroxide  or  (c)  by  the  tricalcium 
aluminate  contained  within  the  hardened  cement  paste.  While  theoretically 
rigorous,  these  expressions  contain  several  uncertainties  which  may  affect  their 
accuracy.  For  example,  the  accessible  porosity  cannot  be  easily  derived.  V/ater 
can  be  stored  in  concrete  in  large  bubbles,  pore*  of  different  radii  (ranging  from 
130  to  500  nm,  see  Barant,  et  al.  1988),  and  intragel  space. The  accessibility  of 
these  passages  depend  on  the  complex  morphology  and  connectivity  of  the  pore 
space  and  can  vary  within  rather  wide  limits  from  case  to  case.  Secondly,  some 
ettringite  is  already  present  (and  some  tricalcium  aluminas  already  consumed)  in 


the  cement  paste  during  hardening  as  a  result  of  the  addition  of  gypsum  needed  to 
avoid  rapid  setting. 

Knowing  the  masses  "formed"  and  "consumed"  during  the  reactions  (4.1) 
and  (4.2),  and  the  densities  of  these  substances,  it  is  rather  straightforward  to 
determine  the  free  volumetric  strain  attributable  to  the  double  decomposition 
reaction.  The  computations  indicate  that  the  mass  of  the  ettringite  formed  will 
almost  never  depend  on  the  calcium  hydroxide  available  in  the  cement  paste. 
Neglecting  that  possibility,  the  volumetric  expansion  computed  from  the 
stoichiometric  (mass  balance)  analyses  can  be  written  as 

ev  =  AV  IV  =  (Mi/ ye)  -  [(Me  /  yc)  +  (Mca  /  ya)\  (4.7) 

From  the  expressions  derived  above  the  free  volumetric  expansion  of  the 
reaction  products  can  be  written  in  a  simple  form  as 

ev  =  3p\  ce  for  ce  <  cs  or  ev  =  3/Jj  cs  =  const  for  ce  ^  cs  (4.8) 


where  Pi  denotes  the  coefficient  of  linear  expansion  of  the  ettringite  crystal. 

The  first  of  the  two  expressions  (4.8)  corresponds  typically  to  the  begin¬ 
ning  of  the  process  during  which  the  reaction  depends  on  the  available  sulphates. 
At  some  point  defined  by 

ce  —  cs  =  (fa  me  Mcm)  /  (<p  ysolV  ma )  ^  ^ 

the  entire  supply  of  the  tricalcium  aluminate  available  in  the  hardened  cement 
paste  is  exhausted  and  the  reaction  is  terminated. 

The  coefficient  P\  is  obtained  comparing  (4.7)  and  (4.8)  in  conjuction  with 

(4.9) 

=  [Mi!  ye)  -  [{Mcc!  yc)  +  {Mca!  ya)\ 

'3fameM cm  /  tyYsolV  ma  (4  JQ) 

However,  the  already  discussed  uncertainties  will  affect  the  accuracy  of  these 
calculations. 
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The  expression  (4. 10)  can  be  further  recast  into  a  more  convenient  form  for 
computations  using  (4.8)2  and  (4.9) 

o  _  ev(ce-  cs)  _  <pySoi  ,me  3mc  ma^ 

1  3 cs  3 me  (ye  yc  ya  (4  10) 

Even  though  the  rigorous  determination  of  the  exact  value  of  the  coefficient 
of  linear  expansion  (3\  (see  also,  Sih,  et  al.  1986)  might  not  fall  within  the  realm 
of  realistic  expectations  in  the  case  of  a  material  with  a  random  morphology  of 
microstructure,  the  fact  remains  that  the  performed  stoichiometric  analysis  clearly 
indicates  the  major  experimentally  determined  trends  in  the  deformation  process. 
In  other  words,  the  expressions  (4.8)  reflect  the  fact  that  the  extent  of  expansion 
is  in  the  short  term  controlled  by  the  concentration  of  the  diffusing  sulphates.  In 
the  long  term  (i.e.  protracted  exposure)  the  volumetric  expansion  depends  only  on 
the  chemical  composition  of  the  cement  (i.e.  the  mass  fraction  of  the  tricalcium 
aluminate  present  in  the  hardened  cement  paste).  The  limit  between  the  short  and 
long  term  exposure  is  defined  by  the  expression  (4.8)  in  conjunction  with  the 
diffusion  equation  (defining  the  time  dependence  of  the  solute  at  a  particular 
location). 


4.2.  Kinetics  of  ettringite  formation 

Even  a  superficial  perusal  of  the  extensive  literature  concerning  the 
formation  of  ettringite  reveals  a  profound  lack  of  consensus  among  the  cement 
chemists  as  to  the  actual  type  of  the  underlying  chemical  reaction  (4.1-2). 
Schwiete,  et  al.  (1966),  Hansen  (1968),  Soroka  (1979)  claim  that  ettringite  is 
formed  by  a  solid-state  (topochemical)  mechanism  at  the  surface  of  anhydrous 
C3A  crystals.  Chatterji  (1968),  and  Mehta  (1976,  1983)  supported  by  scanning 
electron  microscopy  observations  of  the  ettringite  morphology  hypothesized  that 
ettringite  is  formed  by  a  through-solution  mechanism.  Other  authors  (Biczok 
1972,  Ogawa  and  Roy  1982)  Maintain  that  ettringite  may  be  formed  either  by  a 
solid-state  conversion  or  a  through-solution  reaction,  or  both.  If  both 
mechanisms  are  present  the  topochemically  formed  ettringite  is  believed  to  cause 
microcracking.  During  the  crystallization  of  ettringite  from  the  diffusing  solution 


the  crystals  are  likely  to  gradually  fill  (plug)  the  pores  and  initially  even  improve 
the  strength  of  the  cement  paste  (Biczok  1972). 

Naturally,  the  chemical  reaction  irrespective  of  type  is  not  an  instantaneous 
process.  In  the  considered  case  the  situation  is  somewhat  more  complicated  by  the 
fact  that  the  basic  reaction  (4.1-2)  is  heterogeneous,  i.e.  involves  reactants  in 
more  than  one  phase.  Before  the  substances  can  react  they  must  migrate  at  least  to 
the  solid-liquid  interface.  Therefore,  in  addition  to  the  rate  of  chemical  reaction, 
certain  physical  factors  affect  the  rate  of  mass  transfer  between  phases,  and  the 
overall  rate  of  the  heterogeneous  reaction.  According  to  Walas  (1959)  these 
factors  are: 

-  amount  of  the  available  interfacial  surface  to  which  the  rate  may  be  pro¬ 
portional  (surface-controlled  mechanism) 

-  rate  at  which  the  liquid  diffuses  to  and  through  the  interface  film 
(diffusion-controlled  mechanism) 

-  rate  of  the  backward  diffusion  (important  in  reversible  reactions  only) 
Consequently,  a  complete  formulation  of  the  rate  equation  for  a  heterogeneous 
reaction  should  involve  the  rate  equations  governing  both  the  mass  transport  and 
the  chemical  reaction  itself. 

Several  studies  have  addressed  in  the  past  the  kinetics  of  ettringite 
formation.  A  majority  of  these  studies  support  the  view  that  this  process  is 
controlled  by  the  diffusion  of  the  sulphate  ions  to  the  anhydrous  C3A  particles 
(e.g.  Tenoutasse  1968,  Mori  and  Minegishi  1968,  Plowman  and  Cabrera  1984, 
Pommersheim  and  Chang  1986,  1988,  Brown  and  LaCroix  1989). 

Having  in  mind  the  already  mentioned  controversy,  and  in  the  same  time 
trying  to  keep  the  proposed  model  tractable,  it  seems  reasonable  to  introduce  the 
following  simplifications: 

-  diffusion  equations  for  the  sulphate  ions  contained  both  in  MgSC>4  and  in 
newly  formed  CaS04-2H20  will  be  cast  into  a  single  governing  diffusion 
equation.  The  diffusion  will  be  considered  as  a  time  limiting  step  in  the  process 
of  ettringite  formation, 

-  the  through-solution  reaction  (4.1)  is  assumed  to  be  instantaneous  as 
compared  to  (4.2),  and  its  rate  will  not  limit  the  overall  rate, 

-  the  reaction  (4.2)  (considered  to  be  of  the  solid-state  type)  will  affect  the 
overall  rate  but  not  as  much  as  the  diffusion  process.  The  kinetic  equation  for 
(4.2)  will,  nonetheless,  be  included  in  the  governing  equations  of  the  model. 
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In  summary,  it  seems  reasonable  to  consider  the  two  processes 
simultaneously  as  a  diffusion  equation  with  a  sink  (reflecting  "loss"  of  sulphates 
reacting  with  the  active  constituents  of  the  cement  paste).  Moreover,  there  seems 
to  be  no  practical  (or  experimental)  way  to  separate  the  diffusion  of  the  sulphates 
through  the  pores  from  the  diffusion  of  gypsum  through  the  cement  paste  in 
search  for  C3A.  Thus,  these  two  diffusion  processes  are  simply  joined  into  a 
single  one. 

The  diffusion  will  be  considered  in  detail  in  Section  5  of  this  Report.  The 
remainder  of  this  Section  will  be  focused  on  the  formulation  of  the  rate  equation 
for  the  topochemical  reaction  (4.2). 

In  general  the  rate  at  which  the  molar  concentrations  change  during  a 
chemical  reaction  of  the  type  (  pA  +  qB  —>  IL  )  is  defined  by  an  ordinary 
differential  equation  of  the  form  (Walas  1959,  Dawson  1973,  etc.) 

£-*(«-£*)“«> -ft)'  (4.11) 

where:  x  denotes  an  increase  in  the  molar  concentration  of  the  reaction  product 
L,  t  is  time,  a,  b  are  the  initial  molar  concentrations  of  the  reactants  A,  B;  the 
exponents  a  and  /3  are  small  integers  or  fractions.  The  reaction  orders  a  and  p 
are  unity  only  in  the  case  of  simplest  forward  reactions.  In  the  case  of  catalitic, 
crosscatalitic,  autocatalitic  or  reverse  equations  the  order  of  the  reaction  changes. 
The  coefficients  pH  and  q/l  are  determined  by  the  stoichiometry  of  the  reaction. 
In  (4.11)  k  is  the  rate  constant  of  the  process  (depending  on  the  temperature  and 
pressure)  defining  the  fractions  of  all  collisions  that  are  reactive.  Obviously,  the 
dimension  of  this  constant  depends  on  the  overall  reaction  order  (a  +  ft).  It 
should  be  stressed  here  that  there  is  no  generally  useful  theory  that  predicts  the 
rate  laws.  The  rate  law  must  be  found  from  suitably  designed  experiments  (Kotz 
and  Purcell  1987,  Mijovic  and  Ott  1989).  According  to  the  consulted  literature 
such  data  do  not  seem  to  be  available  for  the  considered  double  decomposition 
reaction  (4.1-2).  Therefore,  as  the  simplest  alternative,  it  will  be  assumed  that  the 
topochemical  reaction  (4.2)  proceeds  according  to  a  second-order  rate  law: 

=  *(c°-3c«)(c°-ce) 


(4.12) 
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where  ce  is  the  actual  molar  concentration  of  ettringite,  c£,  are  the  initial 
molar  concentrations  of  gypsum  and  tricalcium  aluminate.  The  equation  (4.12) 
simply  states  that  the  rate  of  ettringite  formation  is  linearly  proportional  to  the 
available  concentrations  of  gypsum  and  tricalcium  aluminate. 

The  rate  equation  (4. 12)  will  be  somewhat  modified  in  the  Chapter  7  where 
molar  concentrations  of  the  substances  will  be  replaced  by  their  mass  percentages 
with  respect  to  some  reference  mass  M.  The  reason  behind  this  conversion  is  the 
solid-liquid  character  of  the  reaction  (4.2).  In  this  case  the  conventional  reference 
volume  (1L  of  solution)  seems  to  be  inappropriate.  Conceptually,  once  the  data 
concerning  the  order  of  the  rate  laws  for  (4.1-2)  become  available  the  proposed 
model  can  be  readily  adjusted. 

Chemical  equilibrium  ("attractor"  state)  is  reached  when  either  all  of  the 
tricalcium  aluminate  is  exhausted  or  when  all  of  the  adsorbed  gypsum  is  spent. 
The  latter  case  will  typically  occur  only  if  the  diffusion  process  is  terminated 
through  the  interruption  in  the  supply  of  the  magnesium  sulphate.  In  both  cases 
the  termination  of  the  reaction  signals  the  end  of  the  expansion  and, 
consequently,  the  cessation  of  the  deformation  process  (Mather  1968). 

It  is  important  to  point  out  that  the  formation  of  an  ettringite  crystal  per  se 
implies  that  it  exerts  large  pressures  on  the  the  surrounding  matrix.  As  argued  by 
Hansen  (1968)  only  some  of  the  crystals  growing  in  places  restricting  their  free 
expansion  (via  constrictions  in  the  pore  size)  will  exert  pressures  on  the  hardened 
cement  paste  sufficient  in  magnitude  to  nucleate  and  subsequently  propagate  a 
crack. 


5.  DIFFUSION  WITH  ADSORPTION 


5.1.  Diffusion  Equation 

As  already  mentioned  the  extent  of  the  considered  chemical  reaction  strongly 
depends  on  the  availability  of  the  magnesium  sulphates. The  magnesium  sulphate 
is  brought  into  contact  with  the  cement  paste  in  form  of  a  solute  contained  in  the 
water  penetrating  through  the  pores  of  the  hardened  cement  paste.  In  general,  the 
conduction  of  water  may  be  ascribed  to  many  causes.  However,  the  pressure  and 
concentration  gradients  are  generally  considered  to  be  the  two  most  frequent 
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generalized  forces  driving  the  water  through  the  solid.  In  this  Report  it  will  be 
assumed  that  the  water  is  driven  through  the  pores  of  the  hardened  cement  paste 
by  the  gradient  of  the  concentration  s(x,  t)  =  grad  c(x,  t)  of  the  solute  in  the  water. 

Even  though  the  derivation  of  the  diffusion-adsorption  equation  is  available 
in  the  literature  a  short  precis  of  the  underlying  theory  seems  to  be  in  order  for 
the  sake  of  completeness  of  this  Report.  Consider  a  problem  of  a  concrete  slab 
saturated  with  inert  water.  At  time  t  =  0  at  least  one  of  the  external  surfaces  of 
the  slab  is  exposed  to  water  with  some  concentration  e  of  a  chemically  aggressive 
solute.  In  the  case  when  the  pressure  gradients  are  small  it  can  be  assumed  that 
the  mixing  (process  of  gradual  homogenization  of  two  fluids)  occurs  via 
diffusion,  i.e.  by  the  molecular  mass  transport  of  solute  from  the  external  fluid  to 
one  contained  within  the  accessible  pores  of  the  slab.  Naturally,  this  process  is 
irreversible  and,  consequently,  energy  must  be  dissipated  during  the  diffusive 
mass  transport. 

The  conservation  of  mass  equation,  in  conjunction  with  the  equation  of 
continuity,  leads  to  the  well  known  partial  differential  equation  governing  the 
process  of  diffusion 

p  ( ^  +  vVc)  =  -div  i  +  q(c,  t) 

dt  (5.1) 

In  (5.1)  c  (x,t )  is  the  concentration  of  the  solute,  v  (x,t)  the  velocity  of 
convective  flow,  p  the  fluid  density  and  q{c,t)  the  source  term  reflecting  the 
change  in  the  concentration  attributable  to  the  adsorption  of  solute  on  the 
adsorbant  (hardened  cement  paste).  The  thermodynamic  current  (flux)  density 
i(c,f)  represents  the  mass  of  solute  diffusing  through  the  unit  surface  during  a 
unit  of  time. 

The  internal  (Gibbs')  energy  density  (Glansdorf  and  Prigogine  1971)  is 
given  by 


de  =  Tds  +  pdV  -  £  p,  dni 
i 


(5.2) 


where  T,  s,  p  and  V  are  the  absolute  temperature,  entropy,  pressure  and 
volume,  respectively.  Additionally,  Pi  and  «,■  are  the  chemical  potential  and 
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the  number  of  particles  (moles)  of  the  i-th  substance  within  volume  V .  From  the 
conservation  of  mass 

X  ni™i  =  1  (5.3) 

where  m,  is  the  mass  of  the  i-th  particle.  Thus, 

Xm‘dni=0  (5.4) 

Expressing,  furthermore,  the  concentration  in  terms  of  the  total  mass  of  particles 
c  =  n\tn\  (5.5) 

the  expression  (5.2)  for  the  internal  energy  density  can  be  rewritten  in  form 


de  -  Tds  +  pdV  -  pdc 

where  (in  the  case  of  mixing  of  only  two  substances) 

^  m  i  m  2 


(5.6) 


(5.7) 


is  the  generalized  chemical  potential. 

The  generalized  chemical  potential  H  is  the  function  of  the  concentration  c, 
temperature  T  and  pressure  p  .  Thus,  its  gradient  can  always  be  written  as 


(5.8) 


where  the  subscripts  to  the  right  of  the  parentheses  denote  variables  held  fixed. 

Assuming  further  the  process  to  be  thermodynamically  stable  and  the  field 
gradients  small,  the  current  density  i  is  (Landau  and  Lifshitz  1988)  a  linear 
function  of  the  gradients  of  the  chemical  potential  and  temperature 

i  =  -aV/i  -  P  VT  (5.9) 
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where  a  and  /3  are  some  constants. 

Substitution  of  (5.7)  into  (5.8)  leads  after  some  obvious  manipulations  to 


i=  ~pD  I  =  -pDyJc 


(5.10) 


where  I  is  the  field  intensity  and 


D  = 


a  l  §££  | 

P  l  dc  )pj 


dp 

dT 


^-(£1 

ap  F  l  dcjp'T 


(5.11) 


the  diffusivity  parameters. 

In  a  general  case  of  viscous  fluids  the  situation  is  much  more  complex  and 
requires  careful  modifications  of  above  expressions  in  the  spirit  of  the  diffusion 
model  suggested  by  Aifantis  (1980).  Restricting  present  discussion  to  non- 
viscous  fluids  the  governing  equation  can  be  obtained  substituting  (5.10)  into 
(5.1) 


+  vVc)  =  div  pd(vc  +  ^VT  +  ^Vp 


+  q{c,T) 


(5.12) 


In  the  present  case  the  temperature  T  is  considered  to  be  a  known  function 
of  time  and  space.  Generally,  this  might  not  be  the  case  and  an  additional  partial 
differential  equation  defining  time  rate  of  change  of  entropy  as  a  function  of 
thermal  flux  (depending  on  T  and  p  )  must  be  added  (Landau  and  Lifshitz 
1988)  to  allow  for  a  simultaneous  determination  of  T  and  c.  At  this  point, 
however,  T  will  be  regarded  as  being  known.  Additional  simplifications  will  be 
introduced  since  in  its  most  general  form  equation  (5.1)  does  not  lend  itself  well 
to  analyses  of  the  type  considered  appropriate  for  the  initial  phase  of  this  research 
project. 

On  the  basis  of  experimental  observations  it  seems  appropriate  to  assume 
that  the  process  is  slow  and  that  the  mass  transport  by  conduction  (i.e.  diffusion) 
vastly  exceeds  that  by  convection.  In  other  words,  it  will  be  assumed  that  the 
second  term  on  the  left  hand  side  of  (5.12)  can  be  neglected. 

Furthermore,  the  solution  concentration  c(x,t)  is,  in  the  considered  case, 
of  the  order  of  several  percent.  Noting  that  the  coefficients  of  thermal  and 
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barrodiffusion  kj  and  kp  must  vanish  when  c  approaches  zero  (since  no 
diffusion  takes  place  at  vanishing  concentration)  the  second  and  third  term  within 
the  brackets  on  the  right  hand  side  of  (5.12)  can  be  neglected  as  well. 

The  source  term  is  typically  represented  in  terms  of  the  time  rate  of  change 
of  the  amount  of  adsorbate  on  the  adsorbent  (hardened  cement  paste)  per  unit 
volume.  At  a  given  temperature  and  pressure  a  definite  relation  is  found  to  exist 
between  the  number  of  molecules  adsorbed  and  the  concentration  of  the  solute. 
This  relation,  defined  by  isotherms,  is  measured  experimentally  for  each  specific 
process  under  controlled  conditions.  When  only  few  sites  on  the  surface  of  the 
pore  are  occupied  the  so-called  Henry's  law  (Compton  1989)  applies  and  the  loss 
of  the  sulphate  concentration  attributable  to  adsorption  is 

cA  =  kic  (5.13) 

where  k\  (kinetic  rate  constant  of  the  process)  is  a  parameter  of  the  system.  In 
some  applications  similar  to  the  presently  considered  one  (Banks  and  Jerasate 
1964)  this  linear  law  was  shown  to  be  quite  useful.  In  other  cases  the  linear  law 
(5.13)  leads  to  a  less  than  satisfactory  fit  of  the  experimental  data  necessitating 
application  of  more  complex  and  nonlinear  isotherms  (Langmuir,  Freundlich, 
etc.).  In  view  of  an  almost  total  lack  of  experimental  data  it  seems  advisable  to 
use  the  simplest  expression  (5.13)  available.  Naturally,  the  expression  (5.13) 
neglects  the  hysteresis  in  the  adsorption-desorption  cycle  which  is,  strictly 
speaking,  appropriate  only  when  the  pores  are  perfect  cylin-ders  and  in  absence 
of  the  residual  liquid  films  covering  the  pore  walls. 

On  the  basis  of  (5. 13)  the  source  term  can  then  be  written  as 

q(c,T,p)  =  -( 1/0 )  ~  =  (5-14) 

With  all  mentioned  simplifications  the  final  form  of  the  diffusion  equation 
(5. 12)  reads 

=  div\pe  (Vc)] 


where 


(5.15) 
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D«  =  D[l  +  <*,/*)]•'  (5.16) 

is  the  effective  diffusivity  accounting  for  the  adsorption  in  the  sense  of 
approximations  (5.13)  and  (5.14).  In  general,  the  diffusivity  is  neither  isotropic 
nor  homogeneous.  These  two  points  will,  however,  be  discussed  in  the  sequel. 

In  the  case  of  constant  diffusivity  and  infinite  half-space  the  Fickian  diffu¬ 
sion  equation  (5. 14)  admits  a  closed  form  solution  in  the  form  of  error  functions 
(Carslaw  and  Jaeger  1959,  Greenkorn  1983,  Landau  and  Lifshitz  1988,  etc.).  For 
a  finite  slab  the  solution  can  be  derived  in  form  of  a  rapidly  converging  infinite 
trigonometric  series  (Carslaw  and  Jaeger  1959).  In  the  case  of  arbitrary  geometry 
and/or  variable  diffusivity  the  solution  can  be  reached  only  through  application  of 
approximate  (such  as  various  types  of  the  method  of  weighted  residuals,  see  for 
example,  Finlayson  1972)  or  purely  numerical  methods.  Since  the  diffusivity  in 
the  considered  case  involves  adsorption  and  microcracking,  neither  of  which  is 
homogeneous,  a  closed  form  analytical  solution  seems  to  be  out  of  reach. 


5.2  Diffusivity 

The  diffusivity  D  is  a  structural  parameter  with  a  strong  dependence  on  the 
tortuosity  of  the  fluid  path  through  the  solid  (Dullien  1979).  In  a  virgin  hardened 
cement  paste  the  diffusivity  depends  on  the  pore  structure,  i.  e.  the  distribution 
(and  connectivity)  of  pores  <P(r)  of  various  widths.  Assuming  the  pores  to  be 
perfectly  random  with  respect  to  their  orientation  and  width  the  virgin 
(undamaged)  hardened  cement  paste  may  be  considered  to  be  of  isotropic 
diffusivity.  For  the  present  purposes  it  will  be  assumed  that  the  accessible 
porosity  <p  and  isotropic  diffusivity  D°  characterizing  diffusion  in  the  virgin 
hardened  cement  paste  are  known. 

As  the  microcracks  within  the  specimen  nucleate  and  grow  the  tortuosity  of 
the  fluid  path  decreases  and  the  diffusivity  increases.  This  problem  was  discussed 
by  many  authors  in  the  past  on  many  levels  of  rigor  and  sophistication 
(Barenblatt,  et  al.  1960,  Salganik  1974,  Dienes  1982,  Englman,  et  al.  1983, 
Wilke,  et  al.  1985,  Bazant,  et  al.  1987,  Babushkin,  et  al.  1987,  Gueguen  and 
Dienes  1989,  etc.). 
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The  purely  continuum  models,  popular  among  analysts  because  of  their 
computational  efficiency,  ignore  the  microscopic  nature  of  the  mixing  process.  In 
general,  these  models  are  based  on  the  volume  averages  of  the  involved  field 
quantities  and  variables  and  are,  therefore,  applicable  in  the  initial  stages  of  the 
process  away  from  the  critical  phenomena  (such  as  percolation  threshold).  Even 
in  these  stages  of  the  process,  i.e.  when  the  microcrack  density  within  the  solid 
is  of  moderate  levels,  it  is  often  argued  that  the  nature  of  the  problem  necessitates 
application  of  statistical  models  (see  Torelli  and  Scheidegger,  1972,  etc.)  based 
on  the  studies  of  a  random  network  of  capillaries.  This  is  especially  important 
when  the  distribution  of  the  capillary  sizes  is  not  isotropic. 

Even  newer  studies  of  the  so-called  percolation  problems  fully  support  the 
stochastic  character  of  the  problem.  Moreover,  these  studies  (see,  for  example, 
Stauffer  1985)  indicate  existence  of  the  scale  effect  (which  is  of  great  interest  in 
the  process  of  interpreting  the  experimental  results)  and  a  radical  change  in  flow 
patterns  and  diffusivity  in  the  neighborhood  of  the  critical  regime  (percolation 
threshold).  These  studies  may,  therefore,  become  pertinent  during  investigations 
of  the  localized  phenomena  characterized  by  the  emergence  of  large  cracks.  This 
work  will,  hopefully,  be  completed  in  the  proposed  continuation  of  this  research 
program. 

The  strategy  adopted  by  Salganik  (1974),  developed  within  the  framework 
of  the  effective  continua  theories,  appears  to  be  most  promising  and  suitable  for 
the  present  purposes.  Consider  a  volume  V  of  an  isotropic  matrix  of  diffusivity 
D°  containing  an  isotropic  inclusion  of  diffusivity  D*.  The  effective  flow  density 
through  the  volume  is  then  given  by  a  vector  equation 


lm  =  (1  +  fjm  (5-17) 

where  f  =  V*  /  V  is  the  volume  fraction  of  the  inclusion.  The  superscripts  (0) 
and  (*)  refer  to  the  matrix  and  the  inclusion.  The  flow  densities  (currents)  are 

jm  ~  DmnSn  jm  =  ^mn  lm  =  (5.18) 


w  here  ^  m  (dc  /  dxm)  is  the  driving  force  (concentration  gradient). 
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Introducing  the  mapping  tensor  Amn  it  is  possible  to  write  the  concentrati¬ 
on  gradient  within  the  inclusion  as  a  function  of  the  effective  concentration 
gradient  in  form 


I  $n  dV  —  A  run  Sm 

Jv1 


(5.19) 


where  V1  is  the  inclusion  volume. 

In  the  present  case  the  interest  is  centered  on  the  determination  of  the  ave¬ 
rage  diffusivity  of  the  hardened  cement  paste  containing  many  spheroidal  in¬ 
clusions.  Moreover,  this  average  diffusivity  is  related  to  some  effective  macro 
continuum  which  in  a  smoothed  manner  relates  to  the  actual  inhomogeneous 
microstructure.  In  order  to  accomplish  the  transition  between  the  effective  macro¬ 
continuum  and  the  inhomogeneous  actual  solid  containing  a  large  number  of 
inclusions  of  different  sizes,  shapes  and  orientations  it  becomes  necessary  to 
introduce  the  concept  of  a  representative  volume  element.  A  representative  volume 
element  must  be  large  enough  to  contain  a  statistically  valid  sample  of 
inhomogeneities.  Selected  in  such  a  manner  it  enables  transiti-on  between  the  two 
scales  mapping  the  actual  solid  on  a  material  point  of  the  effective  continuum. 

For  a  dilute  concentration  of  spheroidal  inclusions  the  expression  (5.19) 
remains  valid  with  a  proviso  that  the  integral  on  the  left  hand  side  of  (5.19)  is 
taken  over  all  inclusions  contained  within  the  representative  volume  element  of 
volume  V  circumscribed  about  the  material  point  in  which  the  diffusivity  is 
sought.  Since  the  number  of  inclusions  within  the  representative  volume  is  by 
definition  large  the  sums  can  be  recast  into  integrals  taken  over  all  sizes  and 
orientations  of  these  inclusions. 

Substituting  (5.18)  and  (5.19)  into  (5.17)  it  becomes  possible  to  derive  the 
expression  for  the  effective  diffusivity  in  form 

Dmn  =  +  {Dmk  -  Dmkj 


Akn  B(Y)  dY 


(5.20) 


Expression  (5.20)  is  subject  to  the  normalizing  condition 
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I 


B(Y)  dY  =  N 


(5.21) 


In  (5.21)  M  is  the  number  of  inclusions  contained  within  the  volume  V.  The 
function  B(Y)  contains  the  probability  density  distributions  of  inclusion  sizes 
and  orientations  within  the  solid  element  dY . 

The  form  of  the  expression  (5.20)  for  the  diffusiv'ty  implies  introduction 
of  the  assumptions  common  to  the  so-called  self-consistent  methods  (ivlura  1982, 
Kunin  1983,  Krajcinovic  and  Sumarac  1987,  etc.).  This  type  of  methods  is  based 
on  the  premise  that  the  direct  interaction  between  inclusions  (and/or  cracks)  is 
negligible  and  that  the  external  fields  of  every  inclusion  do  not  significantly 
differ  from  the  respective  far-fields.  Consequently,  the  val.d.i/  of  the  expression 
(5.20)  is,  strictly  speaking,  limited  to  the  low-to-modest  inclusion  concentration. 
It  is  notable,  that  the  concentration  gradient  s  within  the  inclusion  is 
homogeneous  in  the  case  of  homogeneous  fa’-  field  fluxes  (Landau,  et  al.  1988). 

Assuming  that  both  the  matrix  and  the  inclusion  are  isotropic 

Din  =  -0°  smn  and  =  D*  8mn  (5.22) 

the  relation  between  the  concentration  gradients  in  the  exterior  and  the  interior  of 
the  inclusion  depends  only  on  the  shape  of  the  inclusion  and  can  be  written  in  the 
form  (Shchelchkova  1974,  Salganik  1974,  Landau,  et  al.  1988) 

D°Si  =  D°s ■  +  nik{D'  -  D°)  s*k  (5.23) 

In  (5.23)  riik  is  the  depolarization  tensor  (being  an  analogue  to  the  Eshelby 
tensor  in  solid  mechanics).  The  components  of  the  depolarization  te.isor  can  be, 
in  a  general  case,  determined  only  after  quadratures  of  complicated  integrals. 
Anticipating  application  to  crack-like  defect  it  suffices  to  consider  inclusions  in 
the  form  of  oblate  spheroids,  i.e.  when  lengths  of  the  three  principal  axes  of  the 
spheroidal  inclusion  are 


a  i  =  fl2  =  a 


and 


<23  =  b  =  Pa 


(5.24) 
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In  this  case  the  principal  values  of  the  depolarization  tensor  (see  Landau  and 
Lifshitz  1988)  can  be  derived  in  form  of  simple  expressions 


n(z )  =  1  ± g2  (£  .  tan- ig  )  ,  w(x )  =  n(y )  =  m  =  A-  (1  -  n^z ) ) 

e 3  2 


(5.25) 


where 


e 


1/2 


Therefore,  the  depolarization  tensor  can  be  written  in  form 


(5.26) 


fiik  =  m  Sik  +  (1-  3 m)  8n  8k  3 


(5.27) 


The  expression  (5.22)  generalized  for  all  three  axes  is  then 

j[l-  (1-t]  )m]8ik  -  (  I-77XI  -  3m)<S(-3<5*3  |  =  s-  (5.28) 


or 


(  Aik  )"‘  s*k  =  s[  (5.29) 

where  the  second  rank  tensor  A  is  the  inverse  of  the  bracketed  term  on  the  left 
hand  side  of  (5.28).  The  primes  denote  reference  to  the  local  (principal) 
coordinate  system  of  the  inclusion. 

Thus, 


sk  ^ik  i  (5.  30) 

The  expression  (5.30)  is  written  in  terms  of  the  local  (inclusion)  coordinate 
system  and  must,  therefore,  be  transformed  to  the  global  (specimen)  coordi-nate 
system.  The  general  transformation  rule  for  a  second  order  tensor  is 


Amn  —  &ik  Smi  Snk 


(5.31) 
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where  the  coefficients  of  the  rotation  tensor  g  written  below  in  the  matrix  form 
are  functions  of  the  Euler  angles  <p  and  d 


[5] 


-sin<p  -cosd  cos<p 

cos<(>  -cosd  sin<p 

0  sind 


sind  costp 
sind  sirup 
cosd 


(5.32) 


From  the  expressions  (5.31)  and  (5.32),  in  conjunction  with  (5.28),  after 
some  relatively  simple  but  cumbersome  operations  it  follows  that  the  tensor  A 
acquires  a  relatively  simple  form  of 

Amn  =  A8mn  +  Bnm  nn  (5.33) 


where  the  following  identities  were  utilized 
nm  —  Rmi  and  gm  —  <5,^ 


(5.34) 


The  constants  A  and  B  in  (5.33)  are 


k  A  = 


77  -1 


4  +  (77  -1  )m 


kB  = 


77-I 


1 


1 


77  -  2(77  -  \)m  1  +  (77  -  1  )m 


(5.35) 


In  deriving  (5.35)  the  inclusion  volume  was  taken  as 

V'  =  ±Kaib  =  (5.36) 

with  /J  =  ^  being  the  aspect  ratio  of  the  lengths  of  the  spheroid  axes. 

For  a  low  to  modest  concentration  of  ellipsoidal  inclusions  (such  that  the 
interaction  of  their  external  fields  has  a  negligible  effect  on  the  average  current) 
the  contributions  of  each  inclusion  can  be  merely  superimposed.  In  this  case  from 
(5.20),  (5.22)  and  (5.33)  follows  a  reasonably  simple  expression  for  the  effective 
diffusivity  in  form  of 
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Dik  =  D° 


N 


Sik  +  X  /(')  U(r)  +  B(r)  «i  nk) 
to 


(5.37) 


In  the  case  of  a  large  number  of  inclusions  N  the  addition  in  (5.37)  ceases 
to  be  a  viable  alternative.  Instead,  it  is  necessary  to  rewrite  the  expression  (5.37) 
in  form  of  an  integral 


Dik  =  D° 


ASik  +  Bnirik  )  a3F  sin8d<pddda  db 


(5.38) 


In  (5.38)  the  integration  domain  includes  the  solid  angle 
Y(0<6<n/2  ,  0<<p<2n )  and  the  entire  range  of  the  inclusion  sizes.  In  contrast 
to  (5.20)  expression  (5.38)  implies  that  no  correlation  exists  between  the  sizes, 
aspect  ratios  and  orientations  of  inclusions.  This  seems  to  be  a  reasonable 
consequence  of  the  inherent  randomness  of  the  hardened  cement  paste  structure 
on  the  microscale.  Thus,  the  normalizing  condition  (5.21)  becomes 


J  F(<p,6,a,b)  sind d<p  dO da  db  =  4 n&-  (5.39) 

In  certain  cases  the  sizes  and  orientations  of  the  spheroids  are  not  corre¬ 
lated  allowing  for  a  reasonably  simple  evaluation  of  the  integrals  in  (5.38)  and 
(5.39).  It  is  important  to  notice  that  in  the  present  case  anisotropy  stems  only 
from  the  orientation  of  the  inclusions  since  both  the  inclusions  and  the  matrix  are 
isotropic 

The  evaluation  of  the  effective  diffusivity  (5.38)  is  even  easier  in  the  case 
of  crack-like  defects  for  which  the  parameter  (5  tends  to  zero,  while  t]  =  D*/D° 
tends  to  infinity .  Expanding  the  expressions  (5.25)  into  Taylor  series  it  follows 
that  the  parameter  e  tends  to  1//?.  Thus,  from  (5.35) 


A  =  -B  =  \/p  for  l»p»  1  IT] 


(5.40) 
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which  significantly  simplifies  the  computations. 

In  view  of  the  fact  that  the  inclusions  occupy  only  a  small  part  of  the  entire 
volume  and  in  conjunction  with  the  fact  that  the  local  strain  is  related  to  the 
average  strain  as  the  ratio  of  the  total  volume  to  the  volumes  of  inclusions  it 
seems  reasonable  to  assume  that  the  local  stresses  in  the  vicinity  of  the  expanding 
inclusion  are  much  larger  that  the  average  stresses.  Consequently,  it  follows  that 
the  microcrack  distribution  is  perfectly  random  (with  respect  to  orientations  and 
sizes).  In  the  case  of  isotropic  distribution  of  isotropic  inclusions  embedded  in  an 
isotropic  matrix  the  diffusivity  must  be  isotropic  as  well.  Thus,  the  function  F 
in  (5.38-39)  degenerates  into  a  constant.  In  the  considered  case 

F  =  8  (5.41) 

where  the  bar  above  the  letter  indicate  average  value.  Thus,  after  relatively  simple 
but  laborious  manipulations  the  effective  diffusivity  can  be  derived  from  (5.38)  in 
a  simple  form  of 

Dmn  =  DSmn  (5.42) 


where 


D  =  D°|  1  +22-©) 

with 


(5.43) 


co  = 


N  g3 
v 


(5.44) 


being  the  Budiansky-O'  Connell  (1976)  micromechanical  damage  variable. 

For  the  sake  of  completeness  it  should  be  noted  that  a  different  model  for 
diffusivity  was  suggested  by  Whitaker  (1967)  who  used  the  Taylor  series 
expansion  of  the  product  <v,c>  about  the  vanishing  convective  velocity  <v>  = 
0.  The  result  turns  out  to  be  strikingly,  if  perhaps  fortuitously,  similar  in  form 
since  according  to  the  Whitaker's  model 


Djk  =  D°  (  Sjk  +  Bjk ) 


(5.45) 


where  the  second  rank  tensor  B  represents  the  effect  of  the  tortuosity  (irre¬ 
gularity  of  fluid  paths  and  influence  of  the  cracks,  change  in  pore  radii,  etc.).  An 
even  more  complex  model  using  a  fourth  rank  tensor 
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D  ij  —  dijmn  ‘ 


(5.46) 


explicitly  considers  the  microstructure  of  the  continuum  (through  the  tensor  a) 
and  convective  velocity  v  (Nieman  1969). 

The  proposed  model  defined  by  expressions  (5.38-39)  for  a  general 
distribution  of  cracks,  and  (5.42-44)  for  an  isotropic  distribution  of  cracks,  has  a 
substantial  advantage  in  comparison  to  other  models  which  are  basically 
phenomenological.  Despite  some  uncertainties  with  regard  to  the  determination  of 
the  number  of  cracks  N  and  their  size,  shape  and  orientation  distribution  this 
model  is  firmly  based  on  the  micromechanics  of  the  phenomena.  Therefore,  this 
model  will  readily  allow  future  generalizations  and  extensions  in  quest  of  more 
rigorous  and  sophisticated  analytical  tools.  For  example,  the  model  readily  admits 
considerations  of  tne  influence  of  macrostresses  and  their  sign  (compressive 
versus  tensile)  on  diffusivity  (Salganik  1974),  etc.. 

The  unavoidable  conclusion  of  this  discussion  is  that  a  rigorous, 
deterministic  derivation  of  the  diffusivity  is  not  a  simple  task.  However,  it  is 
undeniable  that  the  influence  of  the  cracks  on  the  diffusion  patterns  is  substantial. 
For  example,  Bazant,  et  al.  (1987)  also  suggest  a  formula  according  to  which  the 
diffusivity  is  proportional  to  the  crack  opening  displacements  and  crack  density 
which  is  in  spirit  identical  to  the  expression  (5.43).  This  formula  was  found  to  fit 
the  experimental  results  rather  well.  Consequently,  it  seems  plausible  to  conclude 
that  the  diffusivity  is,  indeed,  proportional  to  the  microcrack  density,  as 
suggested  by  (5.43).  However,  the  determination  of  the  constant  of 
proportionality  will  at  this  point  be  made  by  fitting  the  experimental  data. 

It  is  finally  important  to  note  that  the  determination  of  the  diffusivity  of 
undamaged  concrete  is  a  problem  in  itself  since  the  diffusivity  of  the  thin 
transition  zones  (at  the  aggregate-cement  matrix  interface)  strongly  influences  the 
flow  patterns.  In  a  virgin  concrete  these  zones  do  not  communicate.  However,  as 
the  cracks  evolve  and  start  propagating  through  the  cement  paste  they  establish 
percolation  paths  which  can  disproportionally  increase  the  macro-diffusivity. 


Thus,  the  effect  of  microcracking  may,  indeed,  be  even  more  pronounced  than 
expected.  This  is  especially  true  in  the  vicinity  of  the  percolation  threshold. 
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6.  STRESS-STRAIN  RELATIONSHIP 


Central  to  the  investigations  of  the  durability  of  concrete  structures  is  the 
determination  of  the  macro-stresses  and  macro-strains  in  the  hardened  cement 
paste  associated  with  the  expansion  of  the  reaction  products  and  the  attendant 
microcrack  nucleation  and  growth  (damage  evolution).  On  the  micro-scale 
chemical  reaction  (4.2)  takes  place  in  small  parts  of  the  volume  in  the 
surrounding  the  pores  with  geometry  conducive  to  the  crystallization  of  ettringite. 
In  general,  the  crystals  will  not  form  in  subcapillary  pores  (in  which  the  surface 
forces  extend  across  the  entire  pore  cross  section,  Hansen  1968).  The  expansion 
pressures  will  not  build  up  in  the  pores  open  to  the  exposed  surface,  either. 
Consequently,  the  crystals  capable  of  nucleating  cracks  will  form  only  in  a  rather 
small  part  of  the  total  volume. 

The  present  task  consists  in  determining  the  compliance  of  the  hardened 
cement  paste  containing  an  ensemble  of  expanding  crystals  and  an  ensemble  of 
microcracks.  At  this  point  it  will  be  assumed  that  all  expanding  crystals  are 
spherical  in  shape  and  that  all  cracks  are  penny-shaped.  Consistent  tc  the 
analyses  in  the  preceding  section,  it  will  be  again  assumed  that  the  defect 
concentration  is  low  to  moderate  rendering  their  direct  interaction 
inconsequential.  It  will  also  be  assumed  that  the  pore  pressures  have  no 
appreciable  effect  on  the  stress-strain  relationship  (see,  for  example,  Biot  1973, 
Rice  and  Cleary  1976,  etc.). 


6.1.  Spherical  inhomogeneities 

Consider  an  isotropic  spherical  inhomogeneity  (occupying  domain  Q  )  with  elastic 
properties  defined  by  a  symmetric  fourth  rank  tensor  C*.  The  free  expansion 
strain  within  the  inhomogeneity,  attributable  to  the  difference  between  the 
volumes  of  reactants  and  the  reaction  product  (ettringite)  is  denoted  by  e'C*  e  Q). 
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The  free  expansion  strain  within  the  inhomogeneity  is,  naturally,  equal  to  one 
third  of  the  volumetric  expansion  (defined  by  the  expression  (4. 8)).  The 

enveloping  matrix  (hardened  cement  paste),  occupying  the  domain  (V'-^2),  is 
also  isotropic,  elastic  and  homogeneous  (in  its  virgin,  undamaged,  state)  with 
material  properties  defined  by  the  symmetric  fourth  rank  tensor  C°.  Consistent 
with  the  previously  employed  notation  the  superscripts  (0)  and  (*)  stand  for  the 
reference  to  the  matrix  and  the  inhomogeneity,  respectively. 

According  to  the  Eshelby's  theorem  the  strains  within  an  ellipsoidal  inho¬ 
mogeneity  are  homogeneous  if  the  far-field  stresses  and  strains  are  homogeneous 
as  well.  Using  this  theorem  and  compensating  for  the  disparity  between  the 
elastic  moduli  tensors  C*  and  C°  of  the  inhomogeneity  and  the  matrix  by  an 
equivalent  uniform  eigenstrain  £*(xe  O)  within  the  inclusion  it  is  relatively 
straightforward  to  determine  the  basic  relations  between  the  stresses  and  strains 
in  the  entire  domain  V .  This  approach  has  been  comprehensively  dealt  with  in  the 
seminal  work  by  Mura  (1982)  rendering  a  recital  of  all  details  within  this  Report 
superfluous.  However,  in  order  to  make  this  Report  self-contained  and  self- 
sufficient  the  major  points  of  the  requisite  analyses,  specific  to  the  problem  at 
hand,  will  be  duly  highlighted  in  the  text  below. 

Denote  by 

e**  =  e*  +  er  (6.1) 

For  an  isotropic,  spheroidal  inhomogeneous  the  eigenstrains  e**  are  related  to  € 
and  far  field  deviatoric  strain  E  °by  (Mura  1982) 

£ij*  =  15(l-v)[cy®  (ji*  -  n)  -  erf  fi*  ]  [(5v-  7)^-  (8  -  10v)/t*  ]~l  (62) 

where  fi  is  the  shear  modulus.  Also,  the  spherical  parts  of  the  strain  tensors  are 
related  by 

eTk  =  3(1  -v)[(K*-K)4k  -  K*eTkk\  [(4v-  2 )K  -  (l+v)/T  ]  l  (6.3) 


where 
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K 


E 

3(1  -  2v) 


(6.4) 


is  the  bulk  modulus  of  the  matrix. 

In  (6.2) 

ey°  =  '  3  S‘i  Qk  (6.5) 

are  the  deviatoric  components  of  the  far-field  strains.  The  moduli  fj.  and  K 
without  superscripts  refer  to  the  matrix. 

The  expression  (6.3)  can  be  further  simplified  if 

{K*  -  K)e°kk  <  K*eTkk  (6.6) 

The  inequality  (6.6)  will  be  satisfied  if: 

-  the  bulk  moduli  of  the  ettringite  and  the  hardened  cement  paste  do  not 
significantly  differ  from  each  other,  and/or  when 

-  the  far-field  strain  is  not  the  result  of  mechanical  loads  but  arises  from 
the  accommodation  of  strains  attributable  to  a  dilute  concentration  of  expanding 
inhomogeneous  inclusions  (beam  or  slab  effect). 

In  absence  of  data  regarding  the  elastic  moduli  of  the  ettringite  it  is  not 
possible  to  estimate  how  much  do  they  actually  differ  from  those  of  the  hardened 
cement  paste.  However,  since  the  dilute  concentration  of  inhomogeneous 
inclusions  is  actually  the  one  examined  within  this  Report  the  local  strains  must 
indeed  be  much  larger  in  magnitude  than  their  averages.  Thus,  from  (6.3)  and 
(6.6)  in  conjunction  with  (4.8)  it  further  follows  for  a  single  ettringite  inclusion 
that 


** 

£kk  ~ 


-3(1  -v]K*  eTkk 
(4v-  2)K  -  (1  +v)K*  f‘ 


-3(1  -vy*  Pi 

(4v-  2)K  -  (1  -t-vyc*  f1 


(6.7) 


where  the  volume  fraction  f1  of  inclusions  is  a  positive  number  not  larger  than 
unity. 


34 


Assuming  further,  for  simplicity,  that  the  elastic  moduli  of  the  inhomoge¬ 
neous  inclusion  (ettringite)  and  the  enveloping  matrix  (hardened  cement  paste)  are 
identical,  it  is  not  difficult  to  deduce  that  in  the  absence  of  the  far-field  stresses 
(from  some  other  source),  the  strains  inside  a  spherical  inclusion  are 

^  =  e9  =  3  V^~v  ^  °e  .  xe  £2  (6-8) 


Outside  the  spherical  inclusion  of  a  radius  a  the  strains  are 


•  -  2  1  +  v  (a) 3  er 
r~  3  1  -  v  Irf  e 


(6.9) 


£e  - 


1  1  4-  V  (a) 3  f,r 

3  1  -  v  'rl  e 


xe(v-fl) 


(6.10) 


where  in  view  of  (4.8) 

=  er  8ij  =  Lev  Sij  (6.11) 

is  the  free  expansion  strain  computed  from  the  stoichiometric  analyses  in  the 
Chapter  4  of  this  Report.  It  is  notable  that  the  inclusion  expansion  can  impart 
large  tensile  hoop  strains  (6.10)  capable  of  crack  nucleation  at  the  inclusion  - 
matrix  interface  (r  =  a  ) 

At  this  point  it  is  tacitly  assumed  that  the  temperatures  generated  during  the 
reaction  are  a  second  order  effect  in  computations  of  the  inclusion  expansion. 
This  point  may  have  to  be  revised  subsequently  since  the  measurements  indicate 
that  a  rather  large  thermal  energy  of  1450J  per  1kg  of  the  reaction  product  is 
liberated  during  a  similar  reaction  taking  place  during  hydration  of  the  cement  in 
the  fresh  cement  paste  (see  Mchedlow  -  Petrosian  1988). 

As  the  solute  diffuses  through  the  specimen  the  chemical  reactions  and  the 
attendant  crystallization  of  the  ettringite  will  take  place  in  a  large  number  of 
locations  surrounding  the  fluid  conducting  pores.  The  exact  sizes,  shapes  and 
locations  of  the  expanding  inclusions  containing  ettringite  crystals  (which  can 
exert  substantial  forces  on  the  surrounding  hardened  cement  paste)  are  generally 
not  known.  Consequently,  it  is  not  possible  to  consider  every  inclusion  and 
superimpose  the  effects  of  all  inclusions  in  the  specimen  to  compute  the  average 
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stress  and  strain  fields.  Similarly  as  in  the  Chapter  5  of  this  Report  this  problem 
is  typically  resolved  within  the  effective  continua  approximation  framework 
discussed  in  great  detail  in  Mura  (1982)  and  Kunin  (1983).  As  was  mentioned 
before  this  approach  has  undeniable  merits  in  the  regimes  preceding  the  critical 
phenomena  but  is,  in  its  conventional  form,  of  limited  utility  thereafter. 

Consider,  therefore,  a  volume  V  containing  a  dilute  concentration  of 
spherical  inclusions  having  a  total  volume  of 


V1  =  f1  V  (6.12) 

The  volume  of  inclusions  is  obviously  proportional  to  the  concentration  of 
ettringite  to  the  third  power 

Since  the  volume  average  of  the  stress  perturbations  (induced  by  inhomo¬ 
geneous  inclusions)  about  their  average  must  by  definition  vanish  the  far-field 
stress  must  be  equal  to  the  average  stress.  It  is  then  possible  to  show  (Mura 
1982)  that  the  average  strain  attributable  to  the  expanding  inclusions  is 


iev)  =/£•*  (6.13) 

The  angular  brackets  in  (6.13)  are  used  to  denote  the  average  (expected)  value  of 
the  bracketed  variable.  Knowing  e'  from  (4.8)  and  (6.11)  it  is  then  possible  to 
compute  from  (6.7)  £**  and  from  (6. 13)  the  average  strain  in  the  material  point. 

Combining  (6.7)  and  (6.13)  it  further  follows  that  the  strain  due  to 
expanding  inhomogeneous  inclusions  is  directly  related  to  the  ettringite 
concentration  ce.  The  spherical  part  of  the  strain  tensor  (6.13)  is  then 

{eick)  =  3  pce  (6.14) 


where  from  (6.7)  and  (6. 13) 


P  = 


-3(1  -  v)K* 

(4v-  2)K  -  (1  +  V)K* 


Pi 


(6.15) 


For  the  case  of  many  inclusions  the  bulk  modulus  K  in  (6.7)  is  that  of  the 
effective  continuum  (matrix  and  the  embedded  inclusions).  Unfortunately,  neither 
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the  self-consistent  scheme  (Budiansky  1965)  nor  the  differential  method  (Hashin 
1988)  do  not  lead  to  an  explicit  expression  for  the  effective  bulk  modulus.  For 
the  present  purposes  it  is,  therefore,  more  convenient  to  apply  either  the  Voigt's 
(upper  bound) 

Kv  =  (l  -/)*0  +  f'K*  (6.16) 

or  Reuss'  (lower  bound) 


_L  =  LzXi  +  JL 
Kr  k°  k* 


(6.17) 


estimate.  In  the  present  case  these  two  estimates  should  not  be  too  far  apart  since 
fl  «  1.  In  (6.16)  and  (6.17)  the  indices  V  and  R  stand  for  Voigt  and  Reuss 
estimates. 

The  Voigt  and  Reuss  expressions  for  the  shear  moduli  are  identical  in  form 
to  ones  for  the  bulk  moduli,  i.e.  (see  Mura  1982) 

Mv  =  (1  -/V  +  /V 

and 

A,  =  LzJL  +  IL 

Vr  /i°  n' 


(6.18) 


(6.19) 


The  elastic  modulus  and  Poisson's 
from  the  expressions 

E  =  2(1  +  v)n 

-  _  3K-2n 
V  ~  2(3 K  +  vl) 


ratio  can  subsequently  be  computed 


(6.20) 


in  either  Voigt’s  or  Reuss'  approximation. 
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6.2  Cracks 

The  second  type  of  defects  consist  of  cracks  of  different  sizes,  shapes  and 
orientations.  For  simplicity  it  will  be  assumed  that  the  cracks  are  all  planar  and 
penny-shaped,  but  of  different  orientation.  In  this  case  a  crack  can  be  modeled  as 
a  oblate  spheroid  of  vanishing  thickness.  The  eigenstrain  in  such  an  inclusion  is 
(Mura  1982,  Krajcinovic  and  Sumarac  1987) 

lint  (a3  4)  =  M3 

a3~*Q  (no  summation)  (6.21) 

In  (6.21)  a3  is  the  vanishing  thickness  of  the  oblate  spheroid  (penny-shaped 
crack),  n  normal  to  the  crack  surface  and  b  the  displacement  discontinuity 
across  the  crack  plane.  Primes  indicate  the  local  (crack)  coordinate  system 
selected  in  such  a  manner  that  the  axis  3  is  collinear  with  the  normal  to  the  crack 
surface. 

Assuming  that  the  crack  during  its  growth  remains  in  the  same  plane  (i.e. 
that  its  growth  is  in  our  case  self-similar)  the  displacement  discontinuity  can  be 
written  in  form  (Hoenig  1978) 

b\  =  [l  -  (x)lajf\l2B'ik  G2kax  (6.22) 

where  a  is  the  crack  radius,  Xj<aj , '  while  B  represents  a  matrix  of  influence 
coefficients  which  depend  on  the  elastic  parameters  of  the  effective  continuum 
(i.e.  original  matrix  and  the  already  existing  cracks).  The  original  matrix, 
naturally,  includes  the  isotropic,  but  not  necessarily  homogeneous,  distributi-on 
of  the  ettringite  inclusions.  Also  o  is  the  external  stress  field  which  is  in  the 
spirit  of  the  self-consistent  method  assumed  to  be  equal  to  the  macro-stresses  a0. 
Rewriting  (6.21)  in  a  symmetric  form  for  the  &-th  crack 


(6.23) 


averaging  the  displacement  discontinuity  over  the  crack  surface 


Abi  = 


(6.24) 


and  introducing  the  coordinate  transformation  relating  the  local  (crack)  and  global 
(specimen)  coordinate  systems 

bt  —  gij  bj  and  O) Ij  —  gimgjn^mn  (6.25) 


(where  the  matrix  [5]  is  defined  by  (5.32))  the  eigenstrain  can  be  shown  to  be  of 
the  following  simple  form  (Krajcinovic  and  Sumarac  1987,  etc.) 


*  (*)  _  n3  h. 
eii  ~  ak 


(jO 

ijmn  umn 


(6.26) 


The  stress  tensor  imposed  symmetries  require 

bijmn  =  Pijmn  +  fiijnm 


(6.27) 


where 


Pijmn  —  ^  Bpq  { 8 pi  8qn  nm  nj  +  8pj  gqn  nm  ni)  (6.28) 

The  transition  to  a  large  number  of  cracks  is  performed  in  the  spirit  of  the 
effective  continuum  models  in  the  manner  discussed  in  the  preceding  chapter  of 
this  Report.  The  total  strain  averaged  over  the  volume  of  the  representative 
element  V  is  then  (Mura  1982) 


(6.29) 


where  the  stresses  and  strains  without  superscript  indicate  perturbations  induced 
by  the  cracks.  Since  the  average  values  of  these  perturbations  vanish,  and  since 
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the  far-field  stresses  and  strains  (denoted  by  the  superscript  ’O')  are  at  the  same 
time  the  expected  values  from  (6.29)  and  (6.26)  it  further  follows  that 

eo  =  [s°  +  £  (aV  V)  b]  :  o°  =  S  :  o°  (6.30) 


where 


S  =  S°  +  £(<z3/V)b  =  S°  +  (^~  bj 


(6.31) 


is  the  effective  compliance.  The  angular  brackets  denote  again  the  volume 
averages.  Assuming  that  no  correlation  exists  between  the  sizes  and  shapes  of  the 
cracks  the  effective  compliance  can  be  written  in  a  simpler  form 


S  =  S°  +  (^)(b(0)) 


(6.32) 


According  to  (6.32)  the  compliance  depends  linearly  on  the  Budiansky-O' Connell 
(1976)  damage  variable  already  defined  by  (5.44).  The  data  related  to  the 
ellipticity  and  orientation  of  the  crack  population  are  stored  in  the  fourth  rank 
tensor  b.  As  in  the  preceding  chapter  in  the  case  of  many  cracks  the  sum  is 
typically  converted  into  an  integral  taken  over  the  probability  density  functions 
defining  the  sizes,  shapes  and  orientations  of  the  cracks  within  the  representative 
volume  elements 


S  =  S° 


I ’«<••*) 


f-B(Y)  dY 


(6.33) 


where  the  normalizing  condition  is  the  same  as  in  (5.21) 


l 


B(Y)  dY  =  N 


(6.34) 


These  integrals  can  in  conjunction  with  some  simplifying  assumptions  be 
computed  in  the  closed  form  (Krajcinovic  and  Fanella  1986,  Sumarac  and 
Krajcinovic  1987,  etc  ).  Assuming  an  isotropic  ensemble  of  penny-shaped  cracks 
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embedded  in  an  isotropic  and  homogeneous  elastic  matrix  the  effective  elastic 
moduli  can  be  derived  in  a  closed,  but  not  explicit,  form  (Budiansky  and 
O'Connell  1976,  Sumarac  1987).  Neglecting  the  influence  of  the  Poisson's  ratio 
the  elastic  modulus  in  a  material  point  of  the  effective  continuum  is 

E  =  (6.35) 

where  a>  is  the  damage  variable  defined  by  (5.44),  while  E  is  the  elastic 
modulus  of  the  material  with  inclusions  but  without  cracks  (6.20). 

The  final  expression  (6.35)  for  the  elastic  modulus  is  valid  for  sufficiently 
dilute  concentrations  of  inclusions  and  microcracks  (i.e.  when  their  direct 
interaction  is  a  second  order  effect). 


6.3.  Stress-Strain  Relationship 

Combining  the  influence  of  the  expanding  ettringite  inclusions  and  already  formed 
cracks  the  relation  between  the  macro-stresses  and  macro-strains  can  be  rewritten 
as 


e°  =  S  :  a0  +  fle**  (6.36) 

The  elastic  modulus  of  the  undamaged  part  of  the  matrix  must  reflect  the  influ¬ 
ence  of  the  ettringite  inclusions  (in  Voigt’s  approximation) 

E°  =  {\  -f')Ecp+flEe  (6.37) 

where  the  subscripts  cp  and  e  refer  to  the  cement  paste  and  ettringite 
crystal. Then, 

o°  =  S1  :{e°-fle”)  (6.38) 

For  the  considered  isotropic  (one-dimensional)  case  the  expression  (6.38)  ac¬ 
quires  a  simple  form 
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°°  -  ®)(£°-^d  (6.39) 

where  the  coefficient  /?  is  obtained  combining  the  expressions  (4. 1 0) i  and  (6. 15). 

The  further  analyses  follow  the  conventional  routes.  The  macrostrains  must 
satisfy  the  equations  of  compatibility  and  requisite  boundary  conditions. 
Additionally,  the  macrostresses  must  satisfy  the  equilibrium  equations  and 
boundary  conditions  on  stresses  as  required  by  the  considered  problem. 


6.4.  Damage  evolution  law. 

The  final  step  in  establishing  the  mechanical  constitutive  theory  is  to  formulate  a 
kinetic  law  relating  the  evolution  of  the  damage  variable  o  as  a  function  of  the 
other  thermodynamic  variables.  Assuming  that  the  cracking  will  occur  in  the 
exterior  of  the  expanding  inclusion  it  seems  reasonable  to  expect  that  the  cracks 
will  grow  out  in  the  direction  of  the  radius  of  the  expanding  ettringite  crystal 
(assuming  it  to  be  spherical  in  shape). 

For  a  two-dimensional  crack-hole  system  shown  in  Fig.  1  an  approximate 
formula  for  the  stress  intensity  factor  reads  (Murakami  1987,  p.264) 


K, 


2 pR 

ix(R+a) 


(6.40) 


where  R  is  the  pore  radius,  a  the  crack  length,  p  the  pressure  within  the  hole. 
In  (6.40)  F0(^)  denotes  a  correction  factor.  For  ( a/R)>  2  it  follows  that  Fo=  1 
and  the  stress  intensity  factor  takes  a  simple  form 


K,  =  2pR 

1  n(R  +  a) 


(6.41) 


For  a  constant  radius  Rq  from  the  equilibrium  in  the  sense  of  Griffith 
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(6.42) 
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where  y  is  the  fracture  (surface)  energy  of  the  hardened  cement  paste.  Assuming 
that  the  fracture  energy  is  not  modified  (reduced)  by  the  adsorption  process  (Rice 
1978)  the  expression  (6.42)  suffices  for  the  determination  of  the  crack  length  a  as 
a  function  of  known  quantities. 


Fig.  1.  Cracking  due  to  expansion  of  an  ettringite  crystal. 


A  conventional  analysis  proves  that 

(unstable)  for  a<R  o  (6. 43)  i 

(stable)  for  a  >  Rq  (6.43)2 

In  other  words,  the  crack  length  vs  pressure  relation  is  discontinuous  a  shown  in 
Fig.  2.  A  very  small  notch  at  the  perimeter  of  the  crystal  will  at  one  point  become 
unstable  and  grow  in  a  single  spurt  until  it  acquires  the  length  a  =  Rq.  After  that 
the  crack  will  increase  its  length  with  the  added  pressure  in  a  monotonic  and 
stable  manner. 
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Fig. 2.  Stability  of  crack  growth 

Since  the  pressure  p  within  the  pore  is  directly  proportional  to  the 
volumetric  expansion  from  (6.42)  it  follows  that  the  crack  length  a  is  proportional 
to  the  square  of  the  ettringite  concentration  ce.  Consequently,  using  (4. 10)  and 
(6.7)  the  damage  parameter  (5.44)  is 

°>  =  zc!  (6  44) 


where 


3  rA  f*.  - - 

/ 1  l'2v  (2  -  4v)  +  (1  +  vfK */K) 


(6.45) 


In  (6.45)  N  denotes  the  number  of  cracks  in  the  reference  volume  V,  v,  E  and  K 
are  the  Poissons’  ratio,  the  Young's  modulus  and  the  bulk  modulus  of  the 
matrix, respectively;  y  is  the  surface  energy  of  the  matrix.  Superscript  refers 
to  the  ettringite. 

It  is  quite  apparent  that  a  very  precise  and  rigorous  determination  of  all  in¬ 
volved  parameters  entering  above  expressions  is  possible  only  in  a  statistical 
sense. Taking  into  the  consideration  that  the  inclusions  are  not  spherical  in  shape 
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and  that  their  size  is  not  deterministic  in  nature  it  is,  at  this  stage,  not  possible  to 
be  more  explicit  then  to  rewrite  the  expression  (6.39)  in  the  form 

0°  =  ^(l  -  g2C$)(e°  -  Pce)  =  E(e°- 1 3  ce)  (6 .46) 

where  gi  and  P  are  parameters  to  be  determined  from  the  expressions  (6.35), 
(6.45),  and  (6.15).  The  expressions  derived  in  this  Report  should  provide  guide 
and  bounds  for  the  determination  of  these  parameters. 


7.  APPLICATION 


The  ultimate  test  of  any  analytical  model  is  in  its  ability  to  replicate 
experimentally  measured  data  and  predict  the  response  of  a  specimen  or  a 
structure.  This  task  is  in  the  present  case  hindered  by:  a  less  than  precise 
documentation  of  the  published  test  results,  unfinished  state  of  the  model  and  the 
physical  reality  (i.e.  essentially  random  nature  of  the  phenomenon).  The 
published  test  results  generally  lack  precise  data  related  to  microstructural 
parameters  (diffusivity,  porosity,  microcrack  density,  etc.)  and  are  rather  vague 
in  defining  what  is  actually  considered  to  be  the  volumetric  strain. 

The  best  set  of  data  related  to  the  expansion  of  mortar  specimens  exposed 
to  a  combined  external  attack  of  MgSOa  and  Na2S04  has  been  provided  recently 
by  Ouyang,  et  al.  (1988)  who  performed  their  experiments  in  concert  with  the 
ASTM  C1012  recommendations.  Mortar  specimens  are  used  in  such  tests  to 
accelerate  the  process  of  evaluating  the  potential  expansion  of  the  concrete 
elements  in  sulphate  solutions  (Price  and  Peterson  1968). 

It  is  commonly  observed  that  in  the  external  sulphate  attack  the  matrix 
cracking  due  to  ettringite  crystallization  accelerates  further  penetration  of  new 
sulphate  ions  which,  in  turn,  intensifies  the  pace  at  which  the  expansive  ettrin¬ 
gite  is  formed.  The  rate  of  the  expansion  is,  thus,  an  increasing  function  of  time, 
as  demonstrated  in  Lea  (1970,  p.352,  Fig.  102)  and  Ouyang,  et  al.  (1988,  Fig.  1). 
Conversely,  in  the  case  of  internal  sulphate  attack  the  amount  of  sulphate  ions 
that  reside  within  the  mortar  decreases  as  the  double  decomposition  reaction  (4. 1- 
2)  proceeds.  Consequently,  the  linear  expansion  versus  time  curves  would  exhibit 
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a  decreasing  slope,  as  reported  by  Ouyang  et  al. (1988).  It  should  be  mentioned 
here  that  the  above  observations  refer  to  mortars  and  concretes  made  of  ordinary 
Portland  cements.  If,  for  example,  silicate  cements  were  used  the  expansion 
curves  in  the  external  sulphate  attack  would  revealed  a  monotonically  decreasing 
slope  (Thorvaldson,  et  al.  1932,  Lea  1970) 

Consider  the  ASTM  C490  test  specimen  used  in  the  measurements  of  the 
length  change,  as  shown  in  Fig.  3a.  Twenty  four  hours  after  casting  the  test  bars 
are  demolded,  then  cured  in  water  while  the  companion  mortar  cubes  are  being 
tested  for  the  compressive  strength.  If  the  required  strength  (20  MPa)  has  been 
attained  the  bars  are  placed  in  a  mixed  solution  of  0. 176  mol/L  MgSCU  and  0. 176 
moI/L  Na2SC>4. 


7 

i 

2.5  cm 

=  25  cm 

- Y 

•s 

£ 

«c 

-- 

JJ^_gage  stud 


Fig.  3  (a)  ASTM  recommended  test  specimen  for  measurements  of  the 
length  change 

(b)  cross  section  of  the  specimen 

The  mix  proportions  of  mortar  and  the  content  of  tricalcium  aluminate  are 
usually  given  as  the  mass  percentages  with  respect  to  cement.  Therefore,  it  is 
necessary  to  convert  the  molar  concentrations  of  sulphates  to  the  mass 
percentages  as  well.  Since  the  analysis  is  confined  to  the  topochemical  reaction 
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(4.2)  the  mass  M  of  V=  1  m3  of  mortar  has  been  selected  as  the  reference1  .  The 
mix  proportions  of  the  mortar  used  in  Ouyang,  et  al.  (1988)  are 

Mcm  •  Msa  :  Mw  =  1:4:  0.6  (7.1) 

where  the  subscripts  cm,  sa,  w  denote  cement,  sand,  and  water,  respectively. 

Assuming  the  average  values  for  the  densities  Jem-  3150  kg/m3, 
jsa  =  2650  kg/m3  (Powers  1968)  the  amounts  of  the  the  constituent  substances  can 
be  computed  from  the  relation 

Mcm+Mja+Mw  _  lm3  (12) 

V  V  V  v  ' 

Jem  Jsa  Jw 

neglecting  the  air  content. 

In  the  considered  case  the  masses  calculated  from  (7.2)  are 

Mcm  =  412  kg  Msa  =  1648  kg  Mw  =  247  kg  ;  M  =  2307  kg  (7.3) 

In  order  to  determine  the  initial  concentration  of  C3A  available  for  the 
reaction  with  the  aggressive  ions  diffusing  through  the  test  specimen  it  is  first 
necessary  to  estimate  how  much  of  the  C3A  anhydride  remains  in  the  cement  after 
the  hydration.  Assuming  that  the  amount  of  gypsum  added  to  portland  cement 
clinker  is  4%  of  Mcm  (Biczok  1972  p.46)  a  true  initial  concentration  of  tricalcium 
aluminate  is 


'•--tH*-  004  fy  <™> 

where:  fa  =  (4.3%,  7%,  8.8%,  12%  )  is  the  C3A  content  in  the  cement  as  given  in 
Ouyang,  et  al.  (1988);  ma  and  mg  are,  as  before,  the  molar  masses  of  tricalcium 
aluminate  and  gypsum.  The  factor  (ma/  3  mg)  results  from  the  stoichiometry  of 
the  reaction  (4. 2). 


1  From  now  on  the  term  concentration  is  to  be  understood  as  the  mass  percentage  unless  stated  otherwise 
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The  next  step  consists  of  computing  the  initial  concentrations  of  sulphates 
with  respect  to  the  total  mass  M .  The  number  of  moles  of  both  MgSC>4  7H2O  and 
Na2SC>4  in  the  saturated  solution  occupying  the  accessible  pores  is 

ns  1  =  ns2  =  ns  =  <p  V  •  0.176  (mol  /  L)  (7.5) 

where  the  subscripts  $1  and  S2  stand  for  MgSC>4-7  H2O  and  Na2SC>4  respectively; 
V=  1  m3,  while  <p  is,  as  before,  the  accessible  porosity. 

The  total  mass  of  sulphates  in  the  considered  volume  becomes 

Ms  =  +  Ms  2  =  n(ms  1  +  mS2)  (7.6) 

where  ms  1  and  msi  are  the  corresponding  molar  masses. 

Taking  in  absence  of  exact  data  0  =  0.4  as  an  average  value  for  mortars 
and  making  use  of  (7.5)  and  (7.6)  the  initial  concentration  of  sulphates  is 


c0  =  =  0.012  (7.7) 

A  tacit  assumption  on  which  the  above  analysis  rests  is  that  the  action  of 
magnesium  sulphate  and  sodium  sulphate  on  mortar  is  similar  in  the  sense  that 
both  substances  serve  as  sources  for  the  aggressive  ions  SO4.  In  practice,  the 
magnesium  sulphate  often  turns  to  be  more  dangerous  because  it  additionally 
attacks  and  decomposes  hydrated  calcium  silicates  (e.g.  Lea  1970).  On  the  other 
hand,  it  has  been  reported  that  a  sodium  sulphate  solution  is  more  corrosive  than 
a  magnesium  sulphate  solution  if  the  calculated  content  of  C3A  in  the  cement 
exceeds  9%  (Biczok  1972  p.50).  In  the  light  of  these  information  the  assumption 
of  an  equal  average  influence  of  both  substances  on  concrete  corrosion  seems 
justified. 

Having  computed  the  initial  concentrations  of  reactants  it  is  now  possible 
to  rewrite  the  rate  equation  for  the  ettringite  formation  (4. 12) 


fTtfi 


-  tA 


d£i 

dt 


(7.8) 
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where  ce  stands  for  the  actual  concentration  of  ettringite,  ac  is  the  initial 
concentration  of  the  adsorbed  gypsum,  a  the  adsorption  constant,  and  c  the  actual 
concentration  of  sulphates. 

A  rigorous  analysis  of  the  double  decomposition  chemical  reaction  (4.1-2) 
would  require  one  kinetic  equation  for  every  substance  entering  the  reaction. 
Expression  (7.8),  and  in  particular  the  assumption  that  the  concentration  of 
adsorbed  gypsum  and  the  sulphate  are  linearly  proportional,  is  an  attempt  to 
simplify  the  analysis. 

The  experimental  measurements  of  the  rate  constant  k  ,  for  the  observed 
chemical  process,  could  not  be  found  in  the  consulted  literature.  Thus,  it  becomes 
essential  to  devise  means  for  a  rational  estimate  of  k  based  on  available 
observations.  For  example,  it  seems  reasonable  to  assume  that  the  rate  at  which 
the  considered  reaction  develops  is  similar  to  the  rate  at  which  the  same  reaction 
proceeds  when  the  water  is  added  to  sand  and  cement  during  the  mixing  of 
mortar.  For  the  latter  case  the  rate  equation  is 


(7.9) 


Knowing  the  initial  concentrations  of  gypsum  added  to  the  clinker  c$  and  the 
tricalcium  aluminate  fa  it  is  then  possible  to  estimate  the  rate  constant  requiring 
that  gypsum  be  consumed  after  a  certain  time  t.  Solving  the  equation  (7.9)  for  k 
we  get 
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(7.10) 


Assuming,  for  illustration,  that  50%  of  eg  will  be  used  up  after  t  =24h  ,  the 
rate  constant  is  computed  from  (7. 10)  to  be  k  =  2.4-10'4  s1.  Naturally,  since  k  = 
^(pressure, T)  the  above  estimate  is  based  on  the  assumption  that  the  temperature 
liberated  by  both  reactions  is  identical  and  that  k  is  not  very  sensitive  to 
pressures  (which  are  much  higher  during  the  formation  of  the  ettringite  in  the 
hardened  cement  paste). 

Finally  it  is  advantageous  to  express  (7.8)  in  a  nondimensional  form  as 
follows 


(7.11) 


=  iciccc  -  0A\2ce){fa  -  0.215  cj 
dt 


where: 


ce  = 


£e 

CO 


H  UkTc° 


(7.12) 


with  Co  given  by  (7.7),  T  =  24  h,  and  fa  defined  by  (7.4). 

Consideration  of  the  diffusion  equation  (5.15)  in  the  actual  three 
dimensional  specimen  requires  costly  and  time  consuming  finite  element 
computations.  Even  though  these  analyses  may  indeed  become  necessary  at  the 
later  stages  of  the  model  development  at  this  point  they  can  hardly  be  justified.  In 
order  to  assess  the  ability  of  the  model  to  replicate  the  experimental  trends  it 
seems  sufficient  to  consider  the  plane  bisecting  the  specimen  in  Fig.  3b  in  halves 
and  assume  the  diffusion  to  be  one-dimensional  (across  the  thickness  only).  In 
this  case  the  diffusion  equation  (5.15)  reduces  to  a  much  more  manageable  form 
of 


dc 

dt 


1+^)| 


where:  co  is  the  damage  variable  (5.44). 

The  initial  and  boundary  conditions  are 

c  (x,0)  =0;  c  (-w0,t)  =  c  (  wo,r)  =  co 
and  additionally  from  symmetry 


(7.13) 


(7.14) 


|(0,<)  =  ° 


(7.15) 


which  assures  a  smooth  curvature  of  the  concentration  profile  near  the  x-axis.  In 
(7.14)  wo  is  the  half-width  of  the  specimen  cross  section.  Making  use  of  the 
relation  (6.44)  the  (7.13)  becomes 
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dc 

dt 


d_ 

dx 


D0(l+glce6)^ 


(7.16) 


The  ettringite  concentration  ce  is  determined  from  the  rate  equation  (7.8).  The 
governing  equations  (7.8)  and  (7.16)  are  highly  nonlinear  and  coupled  rendering 
the  analytical  or  even  approximate  solutions  not  possible  at  this  time.  From  the 
numerical  standpoint  the  analysis  of  the  problem  is  further  complicated  by  the  fact 
that  the  system  of  rate  equations  (7.8)  is  quite  stiff,  to  the  point  that  even  the 
Gear’s  method  (commonly  used  in  this  class  of  problems)  is  not  effective.  At  the 
present  stage  it  seemed  reasonable  to  use  a  somewhat  less  rigorous  (partially 
uncoupled)  approach.  Assume  at  the  very  first  time  step  that  the  diffusivity  D  is 
not  affected  by  the  microcracking  [ci)  =  0  =>  D  =  Z7°).  Thus,  the  solution  of  the 
diffusion  equation  (7.16)  can  be  expressed  in  the  form  of  an  infinite  series 
(Carslaw  and  Jaeger,  1959) 


c{x,f)  =  1  -  -4.  £  2w  Tl  C0S  (  ^  ^  2^ ' ) CXp  -^f~(2n  +  l  J2*  (7.17) 


n  -  0 


where,  consistent  with  (7.12),  the  following  normalization  was  introduced 


D  =  £>°  =  ^ 


(7.18) 


Computing  the  sulphate  concentration  c  from  (7.17)  it  becomes  now  possible  to 
solve  the  rate  equation  (7. 11)  for  the  ettringite  concentration  ce  .  In  the  next  time 
step  this  value  is  inserted  into  the  expression  for  diffusivity  D  =  £>°(l  +  SicM),  a 
new  value  of  c  is  found  from  (7. 17),  and  the  process  repeated.  If  the  time  step  is 
kept  small  this  procedure  is  found  to  lead  to  a  reasonable  solution.  The  numerical 
code  using  the  IMSL  (1987)  routines  for  the  solution  of  a  system  of  21  coupled 
nonlinear  ordinary  differential  equations  (7.11)  at  the  discretization  points  is 
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listed  in  the  Appendix  B.  The  CPU  time  consumed  at  one  run  of  the  program  on 
IBM-3090/XA  was  about  3  minutes. 

The  existing  experimental  data  are  limited  to  measurements  of  the  axial  (or 
volumetric)  strains  averaged  over  the  entire  length  of  the  specimen.  Thus,  the 
estimate  of  the  accuracy  of  the  analytical  model  are  restricted  to  the  comparisons 
of  the  computed  and  measured  average  strains. 

In  the  absence  of  the  externally  applied  tractions  the  average  stresses  in  the 
cross  section  must  vanish 

JfH-0 

a(x,  t)dx  =  0  (719) 

o 


Introducing  further  the  classical  Euler-Bernoulli’ s  hypothesis  of  plane  cross 
sections  the  overall  strain  can  be  computed  from  (6.46)  in  conjuction  with  (7. 19) 


£  =  (3  cQ 


f  ce(x,t)E[  1  -g2(C(fief]dx 
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(7.20) 


where  the  normalized  variables  (7. 12)  and  (7. 18)  have  been  used. 

With  the  gradual  increase  in  the  microcrack  density  the  elastic  modulus 
(6.35)  or  (6.46)  decreases.  At  one  point  the  elastic  modulus  might  vanish,  i.e. 
the  material  becomes  so  heavily  damaged  that  it  looses  its  load  carrying  capacity. 
Consistent  with  the  assumption  that  the  local  stresses  around  expanding 
inclusions  are  dominant  (i.e.  much  larger  than  the  stresses  associated  with  the 
beam  like  behavior)  the  zone  in  which  the  load-carrying  capacity  is  lost  is 
contiguous  with  the  external  surface  (Fig. 4). 

This  type  of  behavior  is  often  modeled  introducing  a  fracture  (or 
percolation)  front,  Fig. 4  (cf.  Kachanov  1986),  separating  the  heavily  damaged 
zone  from  the  material  which  can  still  carry  the  loads.  In  fact  such  a  behavior 
(spalling)  is  commonly  observed  in  the  structural  elements  exposed  to  aggressive 
waters  (e.g.  Biczok,  1972).  The  formulation  of  a  proper  criterion  of  the 
percolation  front  propagation  is  a  complicated  task  considered  to  be  beyond  the 
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scope  of  this  study.  For  the  present  purposes  the  position  of  the  front  will  be 
identified  with  the  plane  in  which  the  elastic  modulus  reaches  the  zero  value. 
From  (6.46)  the  position  of  the  fracture  front  is 

£  =  E[  1  -Si(cf)6]  =  0  (7.21) 

where  cf  is  the  corresponding  critical  concentration  of  ettringite. 


Fig.  4.  Moving  fracture  (percolation)  front. 


Naturally  it  is  well  known  that  the  criterion  E  =  0,  or  (Ocr  =  9/16,  is  but  a 
sign  that  the  accuracy  of  the  self-consistent  model  deteriorates  for  (0  -  0.5  (see 
Hashin  1983).  More  sophisticated  analyses  (Kunin  1982,  Christensen  and  Lo 
1979,  Hashin  1988  )  indicate  that  E  does  not  vanish  at  <y  =  9/16  as  indicated  by 
the  Budiansky-O' Connell  (1976)  analyses.  However,  the  statistical  studies  of 
Hansen,  et  al. ( 1 989)  demonstrate  that  the  failure  does  occur  at  tu  =  0.5  (at  least 
for  the  two-dimensional  lattice  they  considered).  In  other  words,  their  results 
confirm  that  the  crack  interaction,  and  thus,  stress  redistribution  in  the  pre- 
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critical  domain  are  not  crucial  for  the  determination  of  the  critical  point.  Whether 
this  is  actually  true  for  a  three-dimensional  case  cannot  be  assessed  without  a 
serious  study  of  the  corresponding  percolation  process.  At  this  point,  in  view  of 
other  uncertainties  and  the  limited  objective  of  this  Report  it  seems  reasonable  to 
adopt  the  simplest  alternative  and  argue  that  the  criterion  (7.21)  approximates  the 
reality  with  a  satisfactory  accuracy.  The  post-critical  behavior  will  be  reduced  to 
claiming  that  the  stress  in  the  volume  swept  by  the  front  jx|  >  w  vanishes.  In  this 
manner  the  thickness  of  the  narrow  zone  in  which  the  concrete  'softens'  before 
the  stress  actually  vanishes  is  shrunk  to  zero. 

In  view  of  the  condition  (7.21)  the  formula  (7.20)  ceased  to  be  valid  when 
the  percolation  fronts  are  generated.  At  this  point  (7.20)  should  be  replaced  by  an 
approximate  expression  as  follows 
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E(x,t  )ce(  x,t )  dx 
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Jo 


E(x,t)  d x 


(7.22) 


where  ()  stands  for  the  average  value.  In  (7.22)  P  is  the  expansion  coefficient, 
while  w  and  w0  are  specified  in  Fig. 4. 

Results  of  the  performed  computations  are  summarized  in  Figs  5  to  16.  The 
following  numerical  values  were  assumed  for  the  material  parameters  involved  in 
the  model:  D°  =  10‘8  cm2/s  (cf.  Biczok  1972  p .  141),  co  =  0.012  (see  7.7), 

k  =  2.9-10'V1  (7.10),  a  =  0.8,  /?  =  0.1  (4. 10]  and  6.15),  fa  =  { 4.3%,  8.8%,  12%),  v  = 
0.2,  K*IK  =  Q.15.  The  two  fitting  constants  were  assumed  to  be  gi  =  0.51012, 
g2  =  0.71011. 

The  concentration  of  sulphates  across  the  thickness  of  the  specimen  is 
plotted  in  Figs  5  and  6.  For  small  volume  fractions  of  C3A  (or  short  times)  the 
sulphate  concentration  c  is  equal  to  co  only  at  the  external  surface  of  the 
specimen  (x  /  w0  =  1).  However,  if  the  mortar  spalls  off  near  the  surface  the 
heavily  damaged  volume  of  the  mortar  immediately  saturates  as  shown  in  Fig.  6. 

Similar  situation,  Figs  7  and  8,  occurs  with  the  concentration  of  ettringite 
obtained  solving  the  equations  of  the  kinetics  of  the  chemical  reaction  (7. 1 1).  The 
main  difference  between  two  figures  is  in  the  chemical  composition  of  the  mortar 
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(i.e.  mass  of  C3A  ).  The  horizontal  part  of  the  ettringite  concentration  curve  in 
Fig.  8  corresponds  to  the  fractured  zone  in  the  specimen. 

The  change  of  diffusivity,  depicted  in  Fig. 9,  also  shows  the  expected 
patterns.  The  two-fold  increase  of  the  diffusivity  at  the  external  surface  (for  the 
4.3%  C3A  content)  is  of  the  same  order  of  magnitude  as  measured  for  a 
reinforced  concrete  beam  by  Bazant,  et  al.(1987).  For  this  case  of  a  'sulphate 
resistant'  composition  the  change  of  diffusivity  is  restricted  to  the  vicinity  of  the 
external  surface. 

The  reduction  of  the  elastic  modulus,  shown  in  Figs  10  and  11,  develops 
along  the  expected  lines.  The  drop  of  -20%  for  the  mortar  with  12%  C3A  is  in  a 
good  agreement  with  the  experimental  data  from  Ouyang,  et  al.  (1988,  Fig. 3). 
Naturally,  consistent  with  the  adopted  model  the  elastic  modulus  reduces  to  zero 
in  the  part  of  the  volume  swept  by  the  percolation  front  (Fig.  11).  In  the  case  of 
the  sulphate  resistant  mortar.  Fig.  10,  the  damage,  i.e.  drop  in  elastic  modulus  is 
again  confined  to  the  vicinity  of  the  external  surfaces. 

The  stress  distribution  across  the  specimen  thickness  is  plotted  in  Figs  12 
and  13.  For  the  specimen  with  a  safe  mass  fraction  of  C3A  (Fig.  12)  large 
compressive  stresses  occur  near  the  external  surface  while  the  middle  part  of  the 
specimen  is  subjected  to  smaller  tensile  stresses.  As  soon  as  fracture  (spalling) 
occurs  the  stresses  drop  to  zero  and  the  compressed  zone  shifts  inward  (Fig.  13). 
The  computed  values  for  macro-stresses  indicate  that  the  compressive  macro- 
stresses  are  still  below  the  levels  at  which  they  will  influence  the  local  stresses 
around  the  expanding  inclusions.  However,  the  tensile  stresses  in  the  middle 
portion  of  the  specimen  might  indeed  facilitate  cracking  and  accelerate  the  failure. 
This  aspect  of  the  problem  will  be  consideied  in  the  following  stages  of  the 
program. 

The  propagation  of  the  fracture  (percolation)  front  for  the  case  of  a  test  bar 
with  12%  of  C3A  is  plotted  vs  time  in  Fig.  14.  In  absence  of  experimental  data 
this  curve,  while  intuitively  acceptable,  cannot  be  confirmed  in  a  quantitative 
sense.  It  is  notable  that  the  proposed  model  will,  indeed,  allow  for  the 
determination  of  the  safe  duration  of  the  sulphate  attack.  In  the  considered  case 
the  exposures  less  than  60  days  will  not  cause  any  spalling.  Therefore,  once 
thoroughly  tested,  this  model  may  indeed  prove  to  be  a  powerful  design  tool. 
Naturally,  the  model  itself  must  be  improved  by  introduction  of  a  more  realistic 
criterion  for  the  onset  of  critical  phenomena  (percolation  or  spalling). 
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Fig.  15  represents  the  actual  objective  of  this  program,  i.e.  it  tests  the 
ability  of  the  model  in  relating  the  macro-structural,  phenomenological 
observation  (macro-strain  or  expansion)  and  the  chemical  composition  of  the 
mortar.  The  curve  for  the  lowest  concentration  of  tricalcium  aluminate  (4.3%)  is 
duplicated  with  a  remarkable  accuracy.  For  the  mortars  with  a  larger  content  of 
C3Athe  trends  are  qualitatively  duplicated  but  the  fit,  as  expected,  is  somewhat 
less  impressive.  This  is,  naturally,  the  result  of  the  fact  that  the  spalling  which 
occurs  in  these  cases  was  not  adequately  modeled  by  the  present  theory.  A  less 
than  rigorous  selection  of  the  mortar  failure  criterion  (  E  =  0,  i.e.  cu  =  9/16)  (7.21) 
and  the  attendant  uncertainties  with  regard  to  the  post-critical  (softening) 
response  of  the  mortar  resulted  in  a  quantitative  discrepancies  with  respect  to  the 
measured  data.  This  error  is  further  compounded  by  the  fact  that  a  fluid  path 
through  a  system  of  connected  cracks  (percolation  threshold)  will,  in  general, 
occur  well  ahead  of  the  onset  of  fragmentation  (Jaeger,  private  communication 
1990).  Thus,  the  transport  of  sulphates  will  occur  in  a  mode  combining  diffusion 
(through  the  very  thin  pores)  and  percolation.  As  a  result  the  expansion  will  be 
accelerated  establishing  a  closer  fit  between  the  experimental  and  analytical  data. 

In  conclusion,  the  data  plotted  in  Fig.  15  demonstrate  a  remarkable  ability 
of  an  as  yet  incomplete  model  in  replicating  the  trends  in  a  qualitative  way,  and 
provide  a  promise  that  a  better  than  expected  quantitative  fit  can  be  obtained  once 
the  model  is  completed. 

In  a  chemical  sense  the  reactions  (4.1-2)  in  a  mortar  poor  in  tricalcium 
aluminate  consume  all  of  it  without  being  able  to  generate  enough  expansion  to 
endanger  the  integrity  of  the  mortar.  Once  the  C3A  is  totally  consumed  the 
ettringite  crystals  will  cease  to  expand  and  multiply  and  the  mortar  will  remain 
immune  to  further  exposure  to  sulphates  (its  strength  will  be  somewhat  reduced 
by  the  ettringite  inclusions,  see  expressions  (6. 14)  and  (6. 15)  ). 

The  graphs  in  Fig.  16  depict  the  increase  in  the  ettringite  concentration  with 
time  for  mortars  with  three  different  C3A  contents.  It  can  be  seen  that  the 
response  of  mortars  rich  with  C3A  was  dominated  by  the  supply  of  gypsum  (first 
parenthesis  in  7.11).  Thus,  in  the  case  of  continued  supply  of  gypsum  the 
expansion  can  be  continued  until  the  specimen  fractures.  In  contrast  the  response 
of  mortar  with  4.3%  of  C3A  is  limited  by  the  available  tricalcium  aluminate 
(second  parenthesis  in  7.11).  At  one  point,  somewhat  in  excess  of  the  considered 
duration  of  60  months,  the  curve  approaches  a  horizontal  asymptote  and  the 
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specimen  becomes  impervious  to  further  sulphate  attack  (see  Mather  1968, 
Fig.  3). 

The  present  computations  were,  for  the  discussed  reasons,  discontinued 
beyond  £  >  0.15%  since  the  lack  of  information  regarding  the  behavior  of  mortar 
in  the  post  peak  region  (which  occupies  the  ever  increasing  part  of  the  specimen 
volume)  makes  these  computations  incorrect.  At  the  current  state  of  development 
it  was  considered  sufficient  to  demonstrate  the  trends  of  the  response.  Only  the 
analyses  in  the  case  when  the  softening  material  occupies  only  a  small  part  of 
volume  can  be  considered  rigorous.  For  the  case  when  the  softening  material 
occupies  a  large  part  of  the  volume  the  present  model  offers  only  an  approximate 
solution. 

It  should  be  emphasized  again  that  the  model  in  its  present  form  contains 
two  fitting  parameters  g2  and  g \  .In  fact  they  could  be  determined  from  (6.44)  if 
the. involved  micromechanical  constants  are  known. 

The  value  {3  =  0.1  for  the  coefficient  of  linear  expansion  used  in  the 
computations  was  obtained  assuming  (in  absence  of  the  experimental  data)  the 
ratio  of  the  bulk  moduli  of  the  ettringite  and  mortar  to  be  equal  to  0.75.  This 
corresponds  to  the  ratio  of  the  densities  of  two  materials.  It  was  also  accounted 
for  the  fact  that  in  the  analyzed  case  the  solution  contains  both  magnesium  and 
sodium  sulphates  in  equimolar  concentrations.  The  value  for  the  rate  constant  k 
used  in  the  computations  was  lower  than  the  simplified  estimate  discussed  on  the 
page  48.  To  this  end  one  has  to  bear  in  mind  that  the  time  needed  for  the  first  of 
the  double-decomposition  reaction  (4.1-2)  to  complete  has  been  neglected,  while 
the  actual  form  of  the  rate  law  for  the  second  reaction  was  uncertain  as  well. 
However,  the  general  structure  of  the  model  allows  for  the  incorporation  of  any 
modification  once  convincing  and  specific  experimental  data  become  available. 

8.  SUMMARY  AND  CONCLUSIONS 

Despite  the  fact  that  the  theory  is  at  this  point  not  complete  it,  nevertheless, 
strongly  indicates  the  inherent  power  of  microstructural  modeling.  The  major 
advantage  of  this  type  of  modeling  is  in  its  ability  to  relate  the  chemical 
composition  and  microstructure  (porosity,  diffusivity,  etc.)  of  the  cementitious 
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composite  to  its  performance  under  all  conditions.  This  becomes  abundantly  clear 
in  the  case  of  the  investigations  described  and  summarized  in  this  Report. 

As  shown  in  Fig.  15  the  model  clearly  distinguishes  between  the  case  of 
small  concentration  of  tricalcium  aluminate  (which  is  sulphate  resistant  if  less 
than  5.5%,  Ouyang,  et  al.  1988)  and  the  mortars  rich  in  C3A  which  disintegrate 
when  subjected  to  water-borne  sulphate  for  a  prolonged  time  period.  Naturally, 
similar  analyses  can  be  readily  performed  for  an  entire  class  of  related  problems 
such  as  processes  driven  by  pressure  gradients  (flow  through  the  porous  media), 
alcali-aggregate  reactions,  cation-exchange  reactions,  hydrolysis  of  the  cement 
paste,  formation  of  insoluble  and  nonexpansive  calcium  salts,  carbonic  acid 
attack,  frost  action,  etc..  The  apparent  ease  with  which  the  model  succeeds  in 
replicating  the  results  of  a  complex  physico-chemical  process  gives  rise  to 
considerable  confidence  in  its  utility  and  versality,  and  presents  a  strong  impetus 
for  its  future  development. 

As  it  was  duly  noted  in  numerous  places  in  this  Report  the  model  is  by  no 
means  completed.  A  serious  effort  must  be  made  to  clarify  a  number  of 
contentious  arguments  and  remove  some  of  the  uncertainties.  For  example: 

-  The  computations  leading  to  results  plotted  in  Fig.  15  are,  rigorously 
speaking,  correct  only  for  the  curve  with  the  lowest  concentration  of  tricalcium 
aluminate.  Even  though  some  justification  for  the  selection  of  the  fragmentation 
(fracture)  criterion  (to  =  9/16)  is  available  a  much  more  detailed  and  rigorous 
study  of  the  critical  phenomena  and  response  of  the  mortar  in  the  post-critical 
(post-peak  or  softening)  regime  must  be  undertaken  in  order  to  gain  full 
confidence  in  the  model. 

-  The  computation  of  diffusivity  without  considering  the  formation  of  a 
direct  path  through  a  system  of  connected  microcracks  (which  can  happen  at 
rather  low  microcrack  densities  and  much  before  fragmentation)  seriously 
underestimates  the  rate  at  which  the  process  develops  at  the  later  stages  of  the 
deformation.  This  is  quite  obvious  from  the  Fig.  15  for  the  two  mortars  with 
higher  concentration  of  tricalcium  aluminate.  Even  though  the  trends  of 
experimental  data  are  well  replicated  the  model  underestimates  the  expansion  by 
neglecting  the  percolation. 

In  conclusion,  the  proposed  model  formulated  in  the  course  of  this  AFOSR 
sponsored  research  program,  shows  undeniable  promise  and  ability  to  deal  with  a 
complex  and  hitherto  not  well  understood  process.  The  elements  of  the  proposed 
model  are  firmly  based  on  the  chemistry  and  physics  of  the  process,  as  well  as 
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the  properties  of  the  cementitious  composite  and  its  micro-structure.  Once  fully 
developed  and  tested,  this  model  should  contribute  to  the  formulation  of  a  potent 
tool  for  the  optimization  of  the  chemical  composition  of  concrete  in  view  of  its 
purpose  and  the  ambient  to  which  the  structure  will  be  exposed.  However,  in 
order  to  reach  this  stage  additional  basic  analytical  and  experimental  effort  cannot 
be  avoided. 
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APPENDIX  A 


substance 

subscript 

molar  mass 

(g/mol) 

molar  volume 

(cm3/mol) 

density 

(g/cm3) 

Ca(OH)2 

c 

74.0 

33.2 

2.23 

MgS04 

s 

120.4 

72 .1 

1.67 

CaS04-2H20 

g 

18  0 

18.0 

1.00 

H2O 

- 

172.2 

74.1 

2.32 

Mg(OH)2 

- 

58.3 

24.5 

2.38 

3Ca0Al203 

a 

270.2 

88.8 

3.04 

3Ca0Al203-3CaS04-32H20 

e 

1254.5 

725.1 

1.73 

Table  A.l.  Molar  masses  and  molar  volumes  of  reacting  substances 


APPENDIX  B 


77 


:**********«******************************»********»**«***************» 


PROGRAM  TO  SOLVE  A  PARTIALLY  UNCOUPLED  SYSTEM  OF  THE  NONLINEAR 
DIFFUSION  EQUATION  (7.16)  AND  THE  RATE  EQUATION  (7.11)  FOR  THE 
ETTRINGITE  FORMATION,  USING  THE  IMSL  ROUTINE  (DIVPAG). 

SINCE  THE  PROBLEM  IS  STIFF  THE  GEAR'S  BACKWARD  METHOD  IS  ADOPTED 

*****ss*********s*************************************>**********x**** 


IMPLICIT  DOUBLE  PRECISION  (A-H,  O-Z) 

PARAMETER  (NEQ=21,  NPARAM=50,  -NSTEP=50,  LDA=NEQ) 

DIMENSION  A(LDA,NEQ),  CE(NEQ),  PARAM(NPARAM),YM(NEQ),SIGMA(NEQ) 
EXTERNAL  FCN,  FCNJ,  DIVPAG,  DSET,  CONST 
COMMON/PASSE/TEND,XX(NEQ),CC(NEQ),DD(NEQ) 

*******  INITIALIZE  PARAM  -  DIVPAG  ***** 


HINIT  =  2.0 
ENORM  =  1 
METH  =  2 
MXSTEP  =  2000. 

CALL  DSET  (NPARAM,  0.0D0,  PARAM,  l) 

PARAM(l)  =  HINIT 
PARAM(4)=  MXSTEP 

param(io)=  inorm 

PARAM(12)  -  METH 

******  INITIALIZE  OTHER  ARGUMENTS  AND  INTRODUCE  CONSTANTS  ****** 

CALL  DSET(NEQ,  0.0,  CE,  l) 

CALL  DSET(NEQ,0.0,SIGMA,  1 ) 

CALL  DSET(NEQ,0.0,YM,  1 ) 

C  =1.0 
W0  =  1.25 
W  =  W0 
CO  =  0.012 
DO  =  1.0D-8 


******  'SRR'  IS  THE  NORMALIZED  RATE  CONSTANT  'K'  ****** 


SRR  =0.03 
G1  =  03D12 
G2  =  0.7D11 

CECR  =  I.0/(G2**(1.0/6.0))/C0 
CCR  =  1.0 


******  'AL'  IS  THE  ADSORPTION  CONSTANT  'ALFA'  ****** 
******  TA'  IS  DETERMINED  BY  (7.4)  AND  (7.12)  ****** 
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AL  =■  0.8 
FA  =  0.3 
PI  =  CONSTt'PI') 

DELT  A= 1 .0/FLO  AT(NEQ- 1 ) 

BETA  =  0.1 

******  INITIAL  ENTRY  ****** 

IDO  =  1 

******  T  IS  TIME,  TOL  -  TOLERANCE  FOR  ERROR  CONTROL  ***** 


T  =  0.0 
TOL  =  0.1D-3 
TMAX  =  300.0 
DO  10  K=l,  NSTEP 

TEND  =  TMAX*FLO  AT(K)/FLO  AT(N  STEP) 

*****  INTRODUCE  CARSLAW-JAEGER'S  SOLUTION  (7.17)  ****** 

DZ  =  D0*24.0*3600.0/(W**2) 

DO  20  N=1,NEQ 
S  =  0.0 

X  =  1.0-(N-l)*DELTA 
XX(N)=X 


******  COMPUTE  THE  DIFFUSIVITY  ****** 

D  =  DZ*(l.0+Gl*(C0*CE(N))**6) 

DD(N)=D 
DO  30  J=0,500 
R  =  FLOATXJ) 

FI  =((-l)**j)/(2.0*R+1.0)*COS(PI*X*(2.0*R+1.0)/2.0) 

F2  =  EXP(-1.0*D*PI**2/4.0*(2.0*R+1.0)**2*TEND) 

S  =  S+F1*F2 
30  CONTINUE 

C  =  1.0  -  (4.0/Pl)*S 
CC(N)=C 
20  CONTINUE 

******  INTEGRATE  21  ODE  (7.11)  USING  IMSL  ****** 

CALL  DIVPAG  (IDO,  NEQ,  FCN,  FCNJ,  A,  T,  TEND,  TOL,  PARAM,  CE(l)) 
WRITE(8,*)'******************* 

WRITE(8,23)  TEND.T 
23  FORMAT(  1  X,TIME  =  ’,F7.1,2X,F7.l) 

WRITE(8  *)’**«****************' 

WRITE(SJ7) 

77  FORMAT(3X,X,4X,,C,6X,,CE,,7X,D',7X,’W’,6X,'YM’4X;SIGMA'4X, 


******  CHECK  IF  THE  FRACTURE  FRONT  APPEARED  ****** 


DO  66  N=1,NEQ 

IF  (CE(N ).LT.CECR)  THEN 

CONTINUE 

ELSE  IF  (CE(N+1  ).LT.CECR)  THEN 

CE(N)  =  CECR 

CE(N+l)  =  CECR 

CC(N)  =  CCR 

CCXN+1)  =  CCR 

NCR=N 

GOTO  15 

ELSE 

CE(N)  =  CECR 
CC(N)  =  CCR 
ENDIF 

66  CONTINUE 

******  COMPUTE  EXPANSION  (7.20)  OR  (7.22)  USING  SIMPSON’S  RULE 

15  El  =  0.0 
E2  =  0.0 
E3  =  0.0 
E4  =  0.0 

DO  11  N=2,NEQ-lf2 

El  *  El  +  CE(N)*(l.O-G2*(CO*CE(N))**6) 

E3  =  E3  +  l.O-G2*(C0*CE(N))**6 

11  CONTINUE 

DO  12  N=3,NEQ-2,2 

E2  =  E2  +  CE(N)*(  l.O-G2*(C0*CE(N))**6) 

E4  =  E4  +  1.0-G2*(C0*CE(N))**6 

12  CONTINUE 

SMI  =CE(l)*(l.0-G2*(C0*CE(l))**6)+4.0*El+2.0*E2+ 

&CE(21)*(  l.O-G2*(C0*CE(  2 1  ))**6) 

SM2  =  1 .0-G2*(C0*CE(  1  ))**6+4.0*E3+2.0*E4+ 

&1.0-G2*(C0*CE(21))**6 

EX  =  BET  A*C0*((  1 .0- W/ WO  )*CEC  R+W/W0*SM1/SM2) 

******  COMPUTE  THE  YOUNG'S  MODULUS  AND  AXIAL  STRESS  ****** 

DO  100  N=1,NEQ 

YM(N)  =  l.O-G2*(CE(N)*C0)**6 

SIGMA(N)  =  YM(N)*C0*BETA*(SM1/SM2  -  CE(N)) 

100  CONTINUE 

DO  21  N=1,NEQ 

WRITE(8,22)  XX(N),CC(N),CE(N),DD(N),W,YM(N),SIGMA(N)tEX 
22  FORMAT(lX,F4.2,6(F8.5,lX),F9.6) 

21  CONTINUE 


******  COMPUTE  THE  WIDTH  REDUCTION  ****** 


W=  W-NCR*W/FLOAT(NEQ-l) 
10  CONTINUE 


******  RELEASE  THE  WORKSPACE  ****** 

IDO  =  3 

CALL  DIVPAG  (IDO,  NEQ,  FCN,  FCNJ,  A,  T,  TEND,  TOL,  PARAM,  CE(l)) 
END 


****************************************************************** 

****************************************************************** 


******  SUBROUTINE  TO  EVALUATE  THE  RIGHT-HAND  SIDE  OF  (7.11)  ****** 

SUBROUTINE  FCN  (NEQ,  T,  CE,  CEDER) 

IMPLICIT  DOUBLE  PRECISION  (A-H,  O-Z) 

PARAMETER(MEQ=21) 

DIMENSION  CE(MEQ),  CEDER(MEQ) 

COMMON/P  ASSE/TEND,XX(MEQ),CC(MEQ),DD(MEQ) 

SRR  =  0.03 
AL  =  0.8 
FA  =  0.3 
DO  20  1=1, MEQ 

CEDER(I)=  SRR*(AL*CC(I)-0.412*CE(I))*(FA-0.215*CE(I)) 

20  CONTINUE 
RETURN 
END 

*************************************************************** 


******  SUBROUTINE  TO  COMPUTE  THE  JACOBIAN  ****** 

SUBROUTINE  FCNJ(NEQ,  T,  CE,  DCDDC) 

IMPLIQT  DOUBLE  PRECISION  (A-H,  O-Z) 

INTEGER  NEQ 

THIS  SUBROUTINE  IS  NEVER  CALLED  ;  INCLUDED  BECAUSE  THE 
CALLING  SEQUENCE  OF  DIVPAG  REQUIRES  IT. 

RETURN 

END 


